REPORT  DOCUMENTATION  PAGE 


Public  reporting  burden  for  this  collection  of  Information  is  estimated  to  average  1  hour  per  response.  Including  the  time  for  reviewing  Instructions,  searching  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  Information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of 
Information,  including  suggestions  for  reducing  this  burden  to  Washington  Headquarters  Services.  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite 
1204,  Arlington.  VA  22202-4302  and  to  the  Office  of  Management  and  Budget.  Paperwork  Reduction  Project  (0704-01 80),  Washington,  DC  20503. 


1.  AGENCY  USE  ONLY  (Leave  b  2.  REPORT  DATE 

1  -Oct-OO 


4.  TITLE  AND  SUBTITLE 

Modelling  Swell  High  Frequency  Spectra!  and  Wave  Breaking 


6.  AUTHOR(S) 

V.  E.  Zakharov 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESSES) 
Department  of  Mathematics 
617  N.  Santa  Rita  Avenue 

University  of  Arizona 
lucson,  Arizona  1 


3.  REPORT  TYPE  AND  DATES  COVERED 
Annual  Report  1 0/ 1  /99-9/30/00 


5.  FUNDING  NUMBERS 


I N0001 4-98-1 -0070 


8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 


9.  SPONSORING  /MONITORING  AGENCY  NAMES(S)  AND  ADDRESS(ES) 

Office  of  Naval  Research,  Program  Officer  C.L.  Vincent  ONR  321  CD 
Bailston  Centre  Tower  One, 

Arlington,  VA  2221 7-5660 


10.  SPONSORING/MONITORING 
AGENCY  REPORT  NUMBER 


12a.  DISTRIBUTION/AVAILABILITY  STATEMENT 
Approved  for  public  release;  distribution  unlimited. 


12b.  DISTRIBUTION  CODE 


13.  ABSTRACT  (Maximum  200  words) 


In  2000  research  in  framework  of  the  grant  N000  14-98-1  0070  was  performed  by  Dr.  V.  E.  Zakharov  in  collaboration  with  Dr.  A.  Pushkarev 
and  Dr.  A.  Dyachenko.  A  part  of  the  work  was  done  in  collaboration  with  our  French  colleague  Dr.  F.  Dias  and  his  graduate  student  P. 
Guyenne.  The  results  of  the  research  are  summarized  in  five  articles.  Three  of  them  are  accepted  for  publication.  The  results  were 
reported  on  five  international  conferences. 


14.  SUBJECT  TERMS 


15.  NUMBER  OF  PAGES 
3 


1 7.  SECURITY  CLASSIFICATION  1 8.  SECURITY  CLASSIFICATK  T9.  SECURITY  CLASSIFICATION 
OF  REPORT  OF  ABSIRACT 

UNCLASSIFIED  UNCLASSIFIED  UNCLASSIFIED 


16.  PRICE  CODE 


20.  LIMITATION  OF  ABSTRACT 


Standard  Form  298  (Rev.  2-89) 
Prescribed  by  ANSI  Std.Z39-18 


Progress  Report 
N00014-98-1-0070 


"Modelling  Swell  High  Frequency  Spectral  and  Wave  Breaking" 
Principal  Investigator,  V.E.  Zakharov 
1  October  1999  -  30  September  2000 


In  2000  research  in  framework  of  the  grant  N000  14-98-1  0070  was  performed  by  Dr.  V. 

E.  Zakharov  in  collaboration  with  Dr.  A.  Pushkarev  and  Dr.  A.  Dyachenko.  A  part  of  the  work 
was  done  in  collaboration  with  our  French  colleague  Dr.  F.  Dias  and  his  graduate  student  P. 
Guyenne.  The  results  of  the  research  are  summarized  in  five  articles.  Three  of  them  are  accepted 
for  publication.  The  results  were  reported  on  five  international  conferences.  Research  was 
performed  in  the  following  directions: 

1.  Constants  of  motion  in  the  SNL  models. 

The  most  solidly  justified  approach  to  the  description  of  nonlinear  interaction  of  sea 
waves  is  the  use  of  the  kinetic  equation  for  the  spectral  density  of  wave  action  first  derived  by  K. 
Hasselmann  in  1962.  Since  this  time  several  codes  for  the  numerical  solution  of  the  kinetic 
equation  were  developed.  Different  codes  give  qualitatively  similar  but  quantitatively  slightly 
different  results.  So  far  there  was  no  universal  test  making  it  possible  to  estimate  a  quality  of  a 
given  numerical  code.  We  offer  such  a  test  based  on  the  use  of  four  conservative  quantities 
preserved  by  the  kinetic  equation:  wave  action,  energy  and  two  components  of  momentum. 

Using  this  constant  of  motion  is  not  a  trivial  problem.  The  SNL  is  a  complicated 
nonlinear  integral  operator,  describing  the  nonlocal  interaction  of  waves  in  the  K-space.  Due  to 
this  nonlocality  there  is  a  "leakage"  of  the  motion  constant  outside  of  any  finite  domain  of 
integration.  We  made  estimates  of  this  leakage  and  found  that  if  the  domain  of  integration  is 
bounded  by  the  wave  number  K,  the  leakage  is  going  to  the  extended  domain  bounded  by  the 
wave  number  3/2K.  This  result  creates  the  theoretical  foundation  for  development  of  the  "clean" 
universal  test  for  examination  of  different  codes  for  solution  of  the  kinetic  equation.  The  results 
are  summarized  in  the  articles  [1,2]  and  to  be  reported  in  the  Workshop  of  Wave  Prediction  in 
Monterey  (November  6-10  2000). 

2.  Nonlocal  diffusion  models  of  SNL. 

In  1999  we  offered  a  simple  euristic  model  of  SNL  making  it  possible  to  speed  up  its 
computer  simulation  in  four  orders  of  magnitude.  This  model  is  based  on  the  use  of  a  simple 
nonlinear  diffusion  operator.  It  works  very  well  for  broad  in  angle  spectra  but  fails  for  narrow  in 
angle  spectral  distribution.  To  fix  this  point  we  offered  an  array  of  nonlocal  diffusion  models 
including  averaging  the  spectrum  over  angle.  The  new  model  includes  seven  fitting  constants 
which  can  be  found  by  an  optimal  way  as  soon  as  we  will  have  a  "paragon"  exact  numerical 
model  of  the  SNL.  These  results  were  reported  on  the  workshops  of  wave  modeling  in 
Delft(March  2000),  Reykjavik  June  2000)  and  are  to  be  reported  in  Monterey  (November  2000). 
The  article  on  this  subject  is  in  preparation. 
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3.  New  approach  to  the  problem  of  wave  breaking. 

In  collaboration  with  A.  Dyachenko  we  developed  a  completely  new  approach  to  the 
description  of  the  dynamics  of  an  ideal  fluid  with  a  free  surface.  The  new  exact  equation 
describing  the  fluid  dynamics  are  suitable  for  analytical  study  and  numerical  simulatuion.  The 
first  results  of  their  numerical  solution  are  very  promising.  We  hope  by  the  use  of  the  new 
approach  to  soon  develop  a  solidly  justified  analytical  and  numerical  model  of  wave  breaking. 

One  article  is  accepted  for  publication,  another  one  is  in  preparation.  The  result  will  be 
reported  on  the  workshop  on  singularities  in  classical,  Quantum  and  Magnetic  fluids.  20-23 
October,  2000  in  Warwick,  UK. 

4.  One-dimensional  model  of  wave  turbulence  in  deep  water. 

In  collaboration  with  Dr.  F.  Dias  and  P.  Guyenne  we  performed  a  detailed  analytical  and 
numerical  study  of  the  MMT(Maida,  McLaughlin  and  Tabak)  model  of  wave  interaction  in  one 
dimenstion.  We  established  a  fundamental  role  of  the  localized  wave  group(quasisolitons)  in  this 
turbulence.  These  results  can  be  used  for  the  explanation  of  the  formation  of  the  "rogue  waves" 
dangerous  for  coastal  construction. 

They  were  reported  on  the  Conference  of  Wave  Turbulence(June  2000,  Massachusetts) 

Two  articles  are  accepted  for  publication. 
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Turbulence  of  one-dimensional  weakly  nonlinear  dispersive 

waves 

V.  E.  Zakharov,  P.  Guyenne,  A.  N.  Pushkarev,  and  F.  Dias 

ABSTRACT.  The  turbulence  of  weakly  nonlinear  dispersive  waves  is  studied  by 
numerically  integrating  a  three- parameter  one-dimensional  model  equation.  In 
particular  the  validity  of  weak  turbulence  theory  is  assessed.  The  predicted 
power-law  solutions  are  explicitly  determined  and  then  compared  with  the 
numerical  results.  For  both  signs  of  nonlinearity,  it  is  shown  that  the  weakly 
turbulent  regime  is  strongly  influenced  by  the  presence  of  coherent  structures. 

These  are  wave  collapses  and  quasisolitons. 


1.  Introduction 

The  weak  turbulence  theory  developed  by  Zakharov  [8]  is  a  tool  for  obtaining 
the  shape  of  frequency  spectra  in  problems  dealing  with  weakly  nonlinear  dispersive 
waves.  The  applications  of  this  theory  range  from  water  waves  in  hydrodynamics 
to  ion- acoustic  waves  in  plasma  physics.  The  weak  turbulence  theory  is  based  on  a 
hamiltonian  formulation  of  the  problem  where  only  resonant  interactions  between 
weakly  nonlinear  waves  are  taken  into  account.  It  is  then  possible  to  derive  approxi¬ 
mate  equations  by  performing  perturbation  expansions  in  terms  of  the  nonlinearity 
parameter.  Although  the  theory  was  developed  more  than  thirty  years  ago,  few 
proofs,  either  experimental  or  numerical,  have  been  given  to  assess  its  validity  (e.g. 
[7]).  Recently,  Majda  et  al.  [5]  proposed  a  one-dimensional  model  equation  as  a 
basis  to  check  the  validity  of  weak  turbulence  theory.  Numerical  computations  on 
this  model  have  been  reported  in  [1],  [3],  [5]  and  [9].  In  this  paper  we  summarize 
the  most  important  numerical  results  on  this  equation,  which  depends  on  three 
parameters,  and  show  that  the  weakly  turbulent  regime  is  strongly  influenced  by 
the  presence  of  coherent  structures,  namely  wave  collapses  and  quasisolitons. 


1991  Mathematics  Subject  Classification.  Primary  76F55;  Secondary  60H 15. 

Key  words  and  phrases,  weak  turbulence,  Kolmogorov  spectra,  water  waves,  wave  collapse, 
quasisoliton. 
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2.  One- dimensional  model  equation 

The  following  three-parameter  nonlinear  dispersive  equation  was  proposed  by 
Majda  et  al.  [5]: 

(2.1)  i  =  Ukik  +  J  Ti23fc  s(k!  +fo-h-k)  dkidk2dk3 . 

In  equation  (2.1),  which  has  been  written  in  Fourier  space,  denotes  the  k- th 
component  in  the  Fourier  decomposition  of  the  complex  wave  field  and 

(*)  stands  for  complex  conjugation.  Equation  (2.1)  depends  on  three  parameters. 
The  first  parameter,  a,  is  related  to  the  linear  frequency  w*  =  \k\a.  The  second 
parameter,  /?,  is  related  to  the  interaction  coefficient 

(2.2)  Tuak  =  A  | • 

The  third  parameter,  A,  which  also  appears  in  the  interaction  coefficient  (2.2)  and  is 
equal  to  ±1,  governs  the  balance  between  dispersive  and  nonlinear  efFects.  One  can 
use  the  terminology  focusing  for  A  =  —1  and  defocusing  for  A  =  -fl.  The  system 
possesses  two  important  first  integrals,  the  Hamiltonian 

H  =  J u>k\i>k \2dk  +  -  J  Ti23k^i^2^3^k  ^(^i  +  k2  -  k3  -  k)  dkidk2dk3dk 

and  the  wave  action  (or  number  of  particles) 

N  =  J  \M2dk. 

Equation  (2.1)  describes  four- wave  resonant  interactions  satisfying 

(2.3)  ki  +  ki  =  k3  +  k 

(2.4)  4*  ^2  =  u^3  4*  u?k  - 

It  can  be  shown  that  when  a  <  1  the  system  (2.3)-(2.4)  has  nontrivial  solutions 
and  that  dominant  interactions  occur  between  four  waves.  In  all  computations 
the  parameter  a  has  been  set  equal  to  1/2.  This  case  mimics  gravity  waves  in 
deep  water,  whose  dispersion  relation  is  given  by  a;*  =  ( gk )1/2,  where  g  is  the 
acceleration  due  to  gravity.  Computations  for  A  =  +1  were  performed  by  Majda  et 
al.  [5],  Computations  for  A  =  ±1  were  recently  performed  by  Cai  et  al.  [1]  and  by 
Zakharov  et  al.  [9]. 


3.  Kolmogorov- type  spectra 

For  a  weak  nonlinearity,  Zakharov’s  theory  [10]  leads  to  a  kinetic  equation  for 
the  two-point  correlation  function  nk  —  {|^*|2): 

dTXk  f  2 

=  4?r  /  |Ti23*|  {nin2n3  +  nin2nk  -  nin3nk  -  n2n3nk) 

x  S(lo i  4-  ^2  ~  —  Wk)  <fr(&i  4-  k2  —  k3  —  k)  dk\dk2dk3  . 

The  two  main  hypotheses  for  deriving  the  kinetic  equation  are  the  assumptions  of 
gaussianity  and  of  random  phases.  The  stationary  Kolmogorov-type  solutions  are 
given  by 

(3.1) 

(3.2) 


nk 

nk 


ax  |Q|l/34-2/>/3-l+a/3 
aalPI1/3*-2^3-1 
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Table  1 .  Slope  and  flux  sign  for  the  Kolmogorov-type  solutions 
(3.1)- (3. 2).  The  dispersion  parameter  a  is  equal  to  1/2. 


0 

-1 

-3/4 

-1/2 

-1/4 

0 

+3 

power  of  &  in  (3.1) 

-1/6 

-1/3 

-1/2 

-2/3 

-5/6 

-17/6 

sign  of  Q 

+ 

+ 

- 

- 

power  of  k  in  (3.2) 

-1/3 

-1/2 

-1 

sign  of  P 

- 

0 

+ 

+ 

+ 

+ 

and  are  associated  respectively  with  a  particle  flux  Q  and  an  energy  flux  P.  The 
coefficients  oq  and  a<i  denote  the  dimensionless  Kolmogorov  constants.  It  is  impor¬ 
tant  to  emphasize  that  these  solutions  do  not  depend  on  the  sign  of  nonlinearity 
A.  Such  solutions  can  be  written  for  all  values  of  0  and  a  <  1.  But  there  is  a 
physical  argument  which  plays  a  crucial  role  in  deciding  the  realizability  of  the 
Kolmogorov- type  spectra.  Suppose  that  pumping  is  performed  at  some  frequencies 
around  u>fc  =  c Jf  and  damping  at  near  zero  and  uj .  Weak  turbulence  the¬ 
ory  then  states  that  the  energy  is  expected  to  flow  from  uj  to  higher  u^’s  (direct 
cascade  with  P  >  0)  while  the  particles  mainly  head  for  lower  w/c’s  (inverse  cascade 
with  Q  <  0).  Accordingly,  we  need  to  evaluate  the  fluxes  in  order  to  select,  among 
the  rich  family  of  power  laws  (3.1)  and  (3.2),  those  which  are  likely  to  result  from 
numerical  simulations  of  equation  (2.1)  with  damping  and  forcing.  Only  the  cases 

0  <  —3/2  and  0  >  2a  -  3/2 
i.e. 

0  <  —3/2  and  0  >  -1/2  if  a  =  1/2 

are  relevant  because  they  correspond  to  a  particle  flux  towards  large  scales  (Q  <  0) 
and  to  an  energy  flux  towards  small  scales  (P  >  0).  The  signs  of  the  fluxes  are  shown 
in  Table  1  for  a  =  1/2  [9],  Computations  are  performed  in  the  range  0  >  —1/2, 
which  includes  the  case  of  simple  cubic  nonlinearity  (0  =  0)  and  the  case  of  gravity 
waves  (0  —  3) . 


4.  Solitons,  collapses  and  quasisolitons 

The  numerical  results  presented  below  show  that  the  weakly  turbulent  regime 
is  strongly  influenced  by  the  presence  of  coherent  structures.  These  are  solitons, 
quasisolitons  or  collapses.  The  existence  of  solitons  depends  on  the  parameter  A. 
Looking  for  soliton  solutions  of  (2.1)  of  the  form 

4>k(t)  =  e'(n-kV'>t<j>k 
with  and  V  constant  leads  to 

(4.1)  <j>k  =  ~Q  _  fry  +  ~  J  Tuzk  ^1^203  +  &2  -  ka  -  k)  dkidk^dkz . 

For  a  <  1,  the  condition  |&|a  /  0,  Vfc  E  M,  implies  that  the  propagating 

speed  V  is  zero.  Rewriting  equation  (4.1)  in  variational  form: 

6{H  +  QN)  =  0 
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one  can  conclude  that  ‘stationary’  solitons  can  exist  only  if  A  =  —  1.  In  that  case  an 
equilibrium  between  nonlinear  and  dispersive  effects  is  possible.  As  for  nonlinear 
Schrodinger-type  equations,  the  linear  stability  criterion  for  solitons  is  given  by 
dN/dQ  >0  [4].  In  our  case  this  gives 

/3  <  a  -  1 


i.e. 

/?  <  -1/2  if  a  —  1/2  . 

Therefore  the  solitons  are  unstable  in  the  regime  of  interest. 

In  view  of  this  result,  it  is  natural  to  look  at  the  formation  of  collapses.  They 
are  typically  described  by  self-similar  solutions  of  the  form 

where  ^ 

£  =  k  (to  —  p  — - - - ,  e  =  arbitrary  constant . 

2oc 

An  analysis  of  the  convergence  of  the  Hamiltonian  and  of  the  wave  action  integral 
as  t  — y  to  shows  that  necessary  conditions  for  collapses  to  exist  when  a  =  1/2  are 
f3  >  —1/2  for  A  =  —1,  which  coincides  with  the  soliton  instability  criterion,  and 
(3  >  0  for  A  =  +  1.  In  spectral  space,  the  self-similar  solution  behaves  at  t  =  to  like 

(4.2)  nk  ~  k~P+«-2 


which  is  analogous  to  Phillips  spectrum  for  deep  water  gravity  waves  [6] . 

In  the  case  A  =  +1,  quasisolitons  can  exist.  These  are  approximate  solutions  of 
equation  (4.1)  which  look  like  envelope  solitons.  In  the  limit  of  a  narrow  spectrum 
centered  at  k  —  kmi  such  as  Q  —  kmV  +  ^  0,  these  quasisolitons  are  given  by 

t)  ~  <j){x  -  Vi)  e^ikm(x-vt) 

with  Q  and  V  given  by 


m  = 


a(l  —  a)  k 


k P- 

"'771 


a+2 


cosh(/c£) 


,  K  -  —  km  |  <  kjr, 


1 


Q  =  -(1  -  a)  fe"  -  -a  (1  -  a)  k “  2k2  ,  V  =  a  1 . 


When  K/km  is  small,  the  quasisolitons  look  almost  like  true  solitons  and  can  persist 
for  a  long  time.  They  can  play  an  important  role  in  weak  turbulence.  When 
is  large,  the  quasisolitons  can  become  unstable  and  develop  into  wave  collapse. 


5.  Numerical  results 

The  numerical  computations  are  performed  by  adding  to  equation  (2.1)  a  source 
term  in  a  narrow  spectral  band  as  well  as  a  damping  term  containing  a  wave  action 
sink  at  large  scales  and  an  energy  sink  at  small  scales: 

-  ^ .  .  d'ik  -  f  ~  -  - 

(5-1)  l~gf  =  uk^k  +  I  Tu3k  V,i^2^3 £(&i  +  ^2  —  ^3  —  &)  dk1dk2dks 

+  i  {Fk  +  Dk)  i>k 


Fk  =  J2  fjS(k  -  kJ )  and  =  -v-  -  v+  \kf  . 
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Figure  1 .  Level  of  nonlinearity  as  a  function  of  time.  The  param¬ 
eters  are  a  =  1/2,/?  =  0  and  A  =  -f  1  (solid  line);  A  =  —  1  (dashed 
line). 


A  pseudospectral  code  is  used  to  integrate  equation  (5.1).  Details  can  be  found 
in  [9]. 


5.1.  Numerical  results  for  j3  =  0,A  =  ±1.  The  study  is  restricted  to  the 
direct  cascade.  Typical  initial  conditions  are  given  by  random  noise.  Simulations 
are  run  until  a  quasi-steady  regime  is  established  which  is  characterized  by  small 
fluctuations  of  the  energy  and  the  number  of  particles  around  some  mean  value. 
Then  time  averaging  begins  and  continues  for  a  length  of  time  which  significantly 
exceeds  the  characteristic  time  scale  of  the  slowest  harmonic  from  the  inertial  range 
(free  of  the  source  and  the  sink).  In  turn,  the  time-step  of  the  integration  has  to 
provide,  at  least,  accurate  enough  resolution  of  the  fastest  harmonic  in  the  system. 
As  our  experiments  show,  one  has  to  use  an  even  smaller  time-step  than  defined 
by  the  last  condition:  the  presence  of  fast  nonlinear  events  in  the  system  requires 
the  use  of  a  time-step  A t  =  0.005,  which  is  40  times  smaller  than  the  smallest 
linear  frequency  period.  Time  averaging  with  such  a  small  time  step  leads  to  a 
computationally  time-consuming  procedure  despite  the  one-dimensionality  of  the 
problem.  Figure  1  shows  the  time  evolution  of  the  average  nonlinearity  e,  which 
is  defined  as  the  ratio  of  the  nonlinear  part  to  the  linear  part  of  the  Hamiltonian, 
each  part  being  calculated  over  the  whole  field.  Of  course,  this  definition  does  not 
really  make  sense  when  external  forces  are  applied  but  it  provides  a  relatively  good 
estimation  of  the  level  of  nonlinearity  once  the  system  reaches  the  steady  state.  The 
mean  values  of  e  are  0.4  when  A  =  + 1  and  0.2  when  A  =  —  1.  They  are  relatively 
small.  Thus,  the  condition  of  small  nonlinearity  required  by  the  theory  holds  for 
both  systems.  However  the  theory  cannot  explain  the  difference  in  the  values  of  e, 
since  the  same  forcing  is  imposed  in  both  systems. 

The  difference  between  the  focusing  and  the  defocusing  cases  is  even  more 
obvious  when  one  looks  at  the  dissipation  rates  of  particles  and  quadratic  energy 
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Table  2.  a  =  1/2,/?  =  0.  Time-averaged  values  of  the  wave 
action,  quadratic  energy  and  corresponding  fluxes  in  the  stationary 
state. 


A 

N 

E 

Q~ 

Q+ 

P~ 

P+ 

+1 

3 

19 

0.1957 

0.0090 

0.276 

0.258 

-1 

1 

9 

0.0098 

0.0478 

0.014 

1.430 

for  small  wavenumbers: 

Q~  =  2  f  v~\k\~d~  \k?dk,  P~  =  2/  v~\k\-d~  uk\k?dk 

Jk<kf  Jk<kf 

and  for  large  wavenumbers 

Q+  =  2  f  |Ar|d+  |^|2  dk  ,  P+  =  2  f  v*  \k\d+  |^*|2  dk 

Jk>kf  Jk>kf 

where  kf  is  the  characteristic  wavenumber  of  forcing.  Their  time-averaged  values 
in  the  stationary  state  are  collected  in  Table  2. 

The  stationary  isotropic  spectra  of  turbulence  are  displayed  in  Figures  2  and  3. 
Again  the  results  depend  on  the  value  of  A.  For  both  cases  the  theoretical  spectrum 
provides  a  higher  level  of  turbulence  than  the  observed  one.  In  the  focusing  case 
(A  =  —1)  this  difference  is  almost  of  one  order  of  magnitude  but  the  slope  fits 
the  predicted  value  —1  well.  For  A  =  +1,  the  observed  spectrum  almost  coincides 
with  the  weak  turbulence  one  at  low  frequencies  and  then  decays  faster  at  higher 
wavenumbers.  In  this  range,  the  slope  is  close  to  —5/4  as  found  in  [5].  Note  that 
a  new  derivation  of  the  Majda  et  al.’s  spectrum  is  proposed  in  [9]. 

Comparison  of  the  turbulence  levels  and  fluxes  of  particles  Q+  for  both  signs 
of  nonlinearity  leads  to  a  paradoxal  result.  At  A  =  —  1  the  total  number  of  particles 
is  three  times  less  than  at  A  =  H-l,  while  the  dissipation  rate  of  particles  is  higher 
by  one  order  of  magnitude.  It  can  be  explained  only  by  the  presence  in  this  case 
of  a  much  more  powerful  mechanism  of  nonlinear  interactions,  which  provides  very 
fast  wave  particles  transport  to  high  frequencies.  In  our  opinion,  this  mechanism  is 
wave  collapse.  Sporadic  collapsing  events  developing  on  top  of  the  weak  turbulence 
background  could  send  most  of  particles  to  high  wavenumbers  without  violation 
of  energy  conservation,  because  in  each  self-similar  collapse  structure  the  amount 
of  total  energy  is  zero.  Such  a  collapsing  event  is  shown  in  Figure  4.  Note  that 
the  contribution  of  collapses  to  the  high-frequency  spectrum  is  weak  because  they 
produce  a  Phillips-type  spectrum  which  decays  very  fast  as  k  ->  -boo.  In  our  case, 
equation  (4.2)  becomes 

nk  ~  AT3/2 . 

Hence,  only  the  weakly  turbulent  component  k~1  survives  at  large  wavenumbers. 
The  coexistence  of  wave  collapse  and  weak  turbulence  was  also  observed  in  [2]  for 
the  nonlinear  Schrodinger  equation. 
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Figure  2.  /3  =  0,  A  =  —  L  Stationary  and  isotropic  spectra  nk 
vs.  wavenumber.  We  compare  the  computed  spectrum  with  the 
predicted  one  of  Kolmogorov-type  n k  =  ck -1  with  c  =  a2P1^3 
(straight  line). 


FIGURE  3.  j3  =  0,  A  =  +1.  Stationary  and  isotropic  spectra  n* 
vs.  wavenumber.  We  compare  the  computed  spectrum  with  the 
predicted  one  of  Kolmogorov-type  nk  =  ck~l  with  c  =  a2P1^ 
(straight  line). 
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Figure  4.  (3  =  0,A  =  — L  Evolution  towards  collapse  at  x  ~  1 
between  t  =  4999.980  and  t  =  5000.205. 


5.2.  Numerical  results  for  ft  ~  3,  A  =  +1.  At  A  =  -f  1  the  picture  of  turbu¬ 
lence  matches  the  weak  turbulence  prediction  both  quantitatively  and  qualitatively. 
Meanwhile,  the  spectrum  at  high  k's  is  steeper  than  the  theoretical  one.  So  far  we 
cannot  give  a  consistent  explanation  of  this  fact.  We  can  just  guess  that  it  is  some¬ 
how  connected  with  quasisolitons.  As  an  illustration,  Figure  5  shows  the  evidence 
of  the  presence  of  quasisolitons.  More  precisely,  the  system  is  first  separated  into 
several  soliton-like  structures  and  low-amplitude  quasi-linear  waves.  Processes  of 
mutual  interactions  slowly  redistribute  the  number  of  waves  in  a  way  leading  to  the 
growth  of  initially  bigger  quasisolitons  and  the  decay  of  initially  smaller  quasisoli¬ 
tons.  The  final  state  then  consists  of  one  big  quasisoliton  moving  in  a  sea  of  small 
quasi-linear  waves.  This  phenomenon  is  similar  to  the  ‘droplet’  effect  observed  in 
the  non-integrable  nonlinear  Schrodinger  equation  [11].  The  soliton  solution  turns 
out  to  be  a  statistical  attractor  for  the  system:  long  time  evolution  leads  to  the 
condensation  of  the  number  of  particles  into  a  single  soliton  which  minimizes  the 
Hamiltonian. 


6.  Conclusions 

In  conclusion,  the  numerical  results  show  two  types  of  localised  coherent  struc¬ 
tures:  collapses  in  the  focusing  case  (A  =  —  1)  and  quasisolitons  in  the  defocusing 
case  (A  =  -fl).  Their  role  in  the  statistical  properties  of  the  system  can  be  seen 
in  the  spectra.  It  leads  to  a  discrepancy  between  numerics  and  theory.  The  fact 
that  weak  turbulence  is  not  reached  may  be  due  to  the  sparsity  of  resonances  in 
one  dimension  and  to  the  numerical  discretization.  Four- wave  interactions  are  not 
as  efficient,  while  localised  structures  become  dominant.  Therefore  equation  (2.1) 
is  not  such  a  good  model  to  assess  the  validity  of  weak  turbulence  theory. 
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Figure  5.  fi  =  3,  A  =  +1.  Snapshot  of  a  quasisoliton  at  x  ~  3.7 
and  t  —  10880. 
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On  the  dynamics  of  an  ideal  fluid  with  the  free 

surface 

A.  I.  Dyachenko 

Landau  Institute  for  Theoretical  Physics,  2  Kosygin  str.  Moscow,  117334,  Russia 


Abstract 

Exact  cubic  equations  are  derived  for  two-dimensional  flows  of  fluid  with 
the  free  surface.  The  equations  involve  inverse  derivative  of  conformal  mapping 
of  the  domain  occupied  by  fluid  to  the  lower  half-plane,  and  complex  velocity 
of  the  fluid.  They  are  hydrodynamic-type  equations  and  describe  transport  of 
singularities  of  conformal  mapping  in  the  upper  half-plane. 

PACS  codes:  47.10. -fg,  47.15.Hg,  03.40.G,  02.60.Cb 


1  Physical  (natural)  variables 

The  irrotational  motion  of  an  incompressible  inviscid  fluid  with  a  free  surface  is  a 
very  basic  problem  in  fluid  dynamics.  The  equations  describing  the  motion  were  first 
formulated  by  Stokes  in  1845.  In  this  article  two-dimensional  flows  on  the  plane  are 
considered,  that  allows  to  apply  conformal  mapping  theory. 

The  velocity  potential  $  of  the  fluid  satisfies  the  equation 

A(f)  =  0 

in  the  domain  occupied  by  fluid.  The  domain  is  bounded  by  the  free  surface  y  = 
r](x,t)  with  the  boundary  conditions  for  velocity  potential  and  moving  free  surface: 

d<j>  l._,l2 

m  +  5|V01  +OT  =  P> 

on 

+  Vx^x  —  <l>y  (1.1) 

at  y  =  rj(x,  t)  and 


1 


(1.2) 


d$ 

dy 

dj>_ 

dx 


=  0 ,y  ->  -oo, 
=  0,  |a:|  — >  oo. 


Here  g  is  the  gravity  acceleration  and  P  -  constant  pressure  at  the  surface  (Let 
P  =  0).  The  system  (1.1)  and  (1.2)  is  a  Hamiltonian  one  with  the  Hamiltonian 

1  roo  rv(xJ)  1  roo 

%  =  -/  dx  \V(f>\2dy  +  -  [  grj2(x,t)dx  (1.3) 

L  J—  OO  J  — OO  Z  J  —OO 

and  canonically  conjugated  variables 

i>{x,t)  =  (j)(x,t](x,  t),t)  and  rj(x,  t)  (  see  [1]  ): 


drj  _  dH 

dt  dip  ’ 

dip  _  dU 

dt  dg ' 


(1.4) 


Along  with  the  Hamiltonian  there  are  three  more  integrals  of  motion,  amount  of 
fluid: 


d  f°° 

dt  /-oo  ^dx  =  °’ 


and  vertical  and  horizontal  momenta: 


d_ 

dt 

d_ 

dt 


(pydy  —  0, 

(f>xdy  =  0. 


(1.6) 


The  equations  (1.1)  and  (1.2)  are  functionally  nonlinear  and  can  be  hardly  stud¬ 
ied.  The  most  important  known  solutions  were  derived  by  Dirichlet  in  1960,  and 
they  are  ellipsoid  (or  ellipse  in  2D  case),  hyperbola  and  parabola.  Dirichlet  solutions 
are  described  in  detail  in  [2]. 


2  Conformal  variables 

To  simplify  the  equations  (1.1)  and  (1.2)  let  us  perform  conformal  mapping  (see[3|)  of 
the  domain  in  z-plane  occupied  by  fluid  to  the  lower  half-plane  of  a  complex  variable 
w  =  u  +  iv 

— oo  <  u  <  oo,  — oo  <  v  <  0. 
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Now  the  shape  of  the  surface  is  given  parametrically  by  real  and  imaginary  parts  of 
conformal  mapping  given  on  real  axis: 


y  =  y(u,t),  x  =  x(u,t)  =  u  + x(u,t),  (2.7) 

here  x(u,  t)  and  y(u,  t)  are  related  through  Hilbert  Transformation 

y  =  Hx,  x  =  -Hy,  H2  -  -1 


and 

H{f(U)) = p.v.~  r 

7T  J-oo  u'  —  u 

After  conformal  mapping  <j)(x,  y,  t)  ->•  v,  t ).  Let  ^(u,  t)  =  <j){u ,  0,  t),  the  potential 

on  the  surface.  It  was  shown  in[3]  that  y(u,t)  and  ^(u,t)  obey  the  following  system 
of  equations: 

V,  =  (yuH~xu)^  (2.8) 


Hqu  ~ 

+  -j-H^u-gy 


here  g  is  gravity  acceleration,  and 

J  xu  “h  yu  I  +  xu  +  yu. 
Another  form  of  the  equation  (2.8)  and  (2.9)  (see[3])  be 

Vtxu  ~  xtyu  =  -HVU, 

<Htxu-Xt^u  +  gyxu  =  H(^tyu-yt^u  +  gyyu) 


The  integrals  of  motion  acquire  the  following  form: 

'H  =  —  [  ^H^udu  +  -  f  yzxudu, 
J  —  OO  2  J-OO 


9  f°°  _,2„ 


/OO 

yxudu , 

“OO 

/OO 

^xudu, 

-OO 

/OO 

^ yudu . 

-OO 


(2.9) 


(2.10) 


(2.11) 


(2.12) 


When  dealing  with  conformal  mapping  z(w,t)  and  velocity  potential  it  is  conve¬ 
nient  to  write  down  instead  the  real  equations  (2.8)  and  (2.9)  the  complex  equations 
for  z(w,t)  and  complex  velocity  potential  $(u>,f).  One  can  get  these  new  equations 
applying  projector  operator  P 


Hf)  =  j(i + •«)(/) 
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to  the  equations  (2.8)  and  (2.9): 


Zt  iUzu, 

=  iU$u-P&L)  +  ig(z-u).  (2.13) 

\zu\ 

Here  U  is  a  complex  transport  velocity: 


U  =  P{ 


)■ 


It  occurs  that  equations  (2.13)  can  be  simplified  just  by  changing  variables.  In¬ 
deed,  let  us  introduce  instead  of  z(w,t)  and  $(w,t)  another  functions  R(w,t)  and 
V(w,t )  in  a  following  way 


Zyj 


=  -iVzw. 


(2.14) 


(V  is  just  if,  i.e.  complex  velocity).  Note,  that  because  of  z(w,t)  is  conformal 
mapping,  its  derivative  does  exist  in  lower  half-plane  and  does  not  have  zeroes  in 
there.  Thus  function  R(w,  t)  is  analytic  in  the  lower  half-plane  and  has  the  following 
boundary  condition: 


R(w,t) —)■  1,  |it>|  — »  oo,  Im(w)  <  0. 

It  is  obvious  that  boundary  condition  for  V  is: 

y(u;,t)-»0,  |to|  — »  oo,  Im{w)  <  0. 

Then  for  these  analytic  functions  equations  acquire  very  nice  form: 

Rt  =  i(UR'-U'R), 

Vt  =  i(UV'-RP'(VV))+g(R- 1).  (2.15) 

Here 

U  =  P(VR  +  VR). 

The  equations  are  cubically  nonlinear  and  include  linear  integral  operator.  They  are 
similar  to  hydrodynamics  type  equations,  but  given  in  the  complex  plane.  Important 
role  there  plays  complex  transport  velocity  U.  It  conveys  zeroes  of  R(w,t),  which 
are  singular  points  of  conformal  mapping  z(w,t).  Possible  types  of  zeroes  need  to 
be  study  in  the  next  article.  Here  it  should  be  just  mentioned  that  in  the  numerical 
simulations  only  y/w  branch  points  are  observed. 
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It  is  easy  to  include  surface  tension  in  the  equations  (2.15),  one  has  just  to  replace 
the  second  equation  by 


Vt  =  i(UV'  -  RP'(VV)  +  g(R  -  1))  -  2 aRP\Q'Q  -  Q'Q), 

where  a  is  surface  tension  coefficient,  and  Q  =  \/R. 

Equations  (2.15)  keep  the  same  form  for  other  boundary  conditions.  Instead  of 
boundary  conditions  (1.2)  one  can  consider  periodic  boundaries  as  well  as  fluid  of 
the  finite  depth. 

Note  that  equations  include  now  only  derivatives  of  conformal  mapping  and  com¬ 
plex  velocity  potential.  As  regards  integrals  of  motion  they  acquire  more  complicated 
form.  But  if  one  restores  complex  velocity  potential 


then  Hamiltonian  is  equal  to 


n  =  -  Re(<f>)Im(&) 

J—oo 

and  momenta 

/OO  1  roo  1 

Re{$)Re{-)du,  Px  =  /  Re($)Im(-)du. 

-oo  Jth  J -OO  R 

3  Conclusions 

The  main  statement  of  the  article  is  the  following.  The  well-known  equations  for 
flows  of  fluid  with  the  free  surface  (1.1)  are  equivalent  to  the  cubic  equations  (2.15) 
written  for  the  inverse  derivative  of  the  conformal  mapping  R(w,t )  =  l/z’(w,t) 
and  fluid  velocity  V(w,t).  These  equations  are  of  hydrodynamic-type  and  describe 
transport  of  singularities  of  conformal  mapping  in  the  upper  half-plane.  And  moving 
singularities  define  the  shape  of  the  surface. 
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tially  supported  by  INTAS-96-0413  Grant,  by  ONR  Grant  #  N00014-98-1-0070,  by 
RBRF  Grant  No.00-01-00929  and  by  the  Grant  of  Leading  scientific  schools  of  Russia 
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On  conservation  of  the  constants  of  motion  in  the  models  of 
nonlinear  wave  interaction. 

A.Pushkarev,  V.Zakharov 


1.  Introduction 


One  of  the  central  problem  of  the  development  of  the  operational  models  for  sea-wave  prediction  is  an 
adequate  description  of  nonlinear  wave  interaction.  So  far,  the  most  solidly  justified  approach  to  this 
description  is  the  use  of  the  kinetic  equation  for  the  spectral  density  of  wave  action  first  derived  by 
K.Hasselmann  in  1962.  Since  this  time,  several  codes  for  numerical  simulation  of  nonlinear  wave 
interaction  were  developed  (Webb  1978,  Masuda  1980,  Hasselmann  and  Hasselmann  1981,  Resio  and 
Perrie  1991,  Polnikov  1994,  Lavrenov  1998,  Komatsu  and  Masuda  1996,  Van  Vledder  1999,  etc). 

Nonlinear  wave  interaction  is  described  by  a  complicated  nonlinear  integral  operator  and  its 
numerical  simulation  is  a  tricky  problem.  All  existing  algorithms  for  its  simulation  are  cumbersome  and 
time  consuming.  So  far,  they  are  too  slow  to  be  directly  used  in  practical  operational  models  of  wave 
prediction.  Therefore,  the  development  of  faster  approximate  models  of  the  nonlinear  wave  interaction  is  a 
very  urgent  problem. 

The  mostly  common  approximate  model  is  DIA  (Discrete  Interaction  Approximation),  known  also 
as  the  WAM  method.  Hasselmann  and  Hasselmann  offered  it  in  1985.  In  this  model,  the  integral  operator  in 
S  ni  is  replaced  by  a  sum  consisting  of  few  discrete  terms.  Zakharov  and  Pushkarev  proposed  quite  another 
approximate  model,  based  on  the  use  of  the  nonlinear  diffusion  operator  in  1999. 

To  estimate  the  quality  of  an  approximate  model  one  should  compare  its  prediction  with  the  results 
of  numerical  simulation  in  the  framework  of  the  "exact"  kinetic  equation.  To  make  this  comparison  reliable 
on  should  be  sure  that  the  "  exact"  model  is  good  enough  to  be  a  paragon  for  such  test.  Actually,  the  real 
criteria  for  examination  of  quality  of  such  models  are  absent. 

Different  schemes  for  numerical  solution  of  the  Hasselmann  kinetic  equation  give  qualitatively 
similar,  but  quantitatively  slightly  different  results,  and  there  was  no  so  far  a  standard  way  for  estimation  of 
their  reliability.  This  circumstances  makes  the  problem  of  construction  of  fast  approximate  models  of 
nonlinear  wave  interaction  difficult  and  uneasy.  One  cannot  believe  in  an  approximate  model  if  one  cannot 
compare  it  with  a  real  good  standard  for  calibration. 

Meanwhile  there  is  a  natural  way  for  examination  of  numerical  method  for  solution  of  the  wave 
kinetic  equation.  This  is  the  control  of  conservation  of  the  basic  constants  of  motion  -  wave  action,  energy, 
and  two  components  of  momentum.  Similar  approach  is  widely  used  in  applied  mathematics  in  the  case 
when  physical  situation  is  described  by  conservative  ordinary  or  partial  differential  equations. 

In  the  case  of  nonlinear  wave  interaction,  the  situation  is  more  complicated.  It  is  described  by  not 
differential,  but  by  the  integral  operator,  which  is  non-local  in  the  k-space.  Any  scheme  of  numerical 
integration  of  the  kinetic  equation  operates  in  some  finite  domain  of  this  space,  always  bounded  in 
frequency  and  sometimes  limited  in  angle.  Integrals  of  motion,  contained  in  any  bounded  domain  are  not 
conserved;  the  non-linear  wave  interaction  carries  them  out  of  the  domain.  Due  to  non-locality  of  the 
Sni  operator,  this  leakage  cannot  be  interpreted  just  as  a  flux  through  a  boundary  of  the  domain.  The  loss 

of  the  motion  constants  from  finite  domains  is  not  a  mathematical  abstraction.  Its  is  a  real  and  a  very  strong 
physical  effect.  In  many  cases,  the  transport  of  motion  constants  is  the  major  mechanism  defining  shape  of 

spectra.  For  instance,  a  typical  asymptotic  behavior  of  energy  spectrum  at  large  frequency  e,  -  r4  ^ 

the  result  of  constant  transport  of  wave  energy  to  the  large  frequency  region. 

In  this  article,  we  propose  a  modified  method  of  calculation  of  conservation  of  the  motion  constant 
in  finite  domain  making  possible  to  take  into  consideration  the  leakage  outside  a  domain.  We  call  it  "clean 


test",  which  allows  accurate  estimation  of  the  quality  of  any  algorithms  for  numerical  solution  of  the  wave 
kinetic  equation. 


2.  Are  the  constants  of  motion  really  conserved? 


In  the  absence  of  pumping  from  the  wind  and  dissipation  the  Hasselmann  equation  reads: 


dn  „ 

dl~bnl 

(2.1) 

$ni  =  Jl ^kkik2k}  ^  ^(k  +  k{  —  k2  —  k3)S(cok  +  0)ki  +  0)h  +  COk} ) X 
(\nk2nk}  +nknk2nk}  -nknknk2  -nknknki)dkxdk2dk3 

(2.2) 

is  considered  that  equation  preserves  the  following  constants  of  motion 

N  =  J  ndk  -  Wave  action 

(2.3) 

E  =  J  (Qknkdk  -  Energy 

(2.4) 

M  =  J  knkdk  -  Momentum 

(2.5) 

Are  these  constants  really  constant?  To  prove  conservation  of  these  integrals,  one  must  prove  validity  of 
following  identities: 


dN 

dt 

)Snldk=  o 

(2.6) 

dE  t 

irJ 

cokSnldk  =  0 

(2.7) 

dM 

r  r  — 

dt  ~ 

J  kSnldk  =  0 

(2.8) 

These  identities  are  trivial  if  one  can  change  the  order  of  integration  by  different  kt .  If  this  operation  is 
possible,  one  can  transform,  for  instance,  the  expression  (2.7)  to  the  form 

J MkSnA  =  f (G)k+  (0K  +  cok  +  CO.  )S(k  +  kl-k2-k3) 

cv  (2.9) 

S(0)k  +(Qki+Q)ki  +  )nknknkdkxdk2dk3 

As  it  is  known  from  the  classical  calculus,  the  operator  of  change  of  the  integration  order  in  improper 
integrals  is  allowed  if  the  integrand  decays  fast  enough  at  infinity.  Let  us  consider  this  question  in  detail. 

In  equation  (2.2),  as  well  as  in  formulae  (2.3)-(2.9),  the  integration  is  going  on  the  infinite  domain. 
In  reality,  both  in  experiment  and  in  computer  modeling  the  domain  of  integration  is  finite.  Thus  to  check 
the  identities  (2.6)-2.8)  we  should  first  consider  a  finite  domain.  This  is  a  quite  nontrivial  procedure. 
Suppose  that  the  domain  of  integration  is  finite 


\k\<p 


(2.10) 


One  can  denote 


n 


p 

k 


nk,\Jc\<  p 
0, 1  k  l>  p 


(2.11) 


By  plugging  Dk  instead  of  Dk  into  S  nt  one  get  by  definition  S n>  — >  .  Apparently  integrand  in  Sk° 

has  bounded  support  and  change  of  order  of  integration  is  permitted  at  any  value  of  p.  Hence 

js^dk=  jslf>dk+  js<fdk  =  0 

\k\<p  U:l>/? 

jotSl!>dk=  j<o,Sl!'dk  =  0  (2  12) 

Ifckp  l£l>p 

f  kS(fdk  =  [  k(okS[p)dk  +  f  kS(fdk  =  0 

\k\<p  \k\>p 


Let  us  denote 


Np  =  J ndk  ,  Ep  =  jo)knkdk,  Mp  =  jknkdk 


\k\<p  \k\<p  \k\<p 

Now  one  can  find  balance  of  the  motion  constants  in  the  domain  | k  \c  /Z 

4  jn,dk=  j S^dk  -  -  | S'fdk  =  -Q(p) 


dt  dt 


\k\<p 


\k\<p 


M>p 


3E""  =  4  J =  j °>X!'dk  =  -  | <atS<fdk  =  -P(p) 


dt  dt 


\k\<p 


l&lcp 


\k\>p 


\kn*dk  =  \kS^dk=-\kSfdk=-K(p) 


dt  dt 


\k\<p 


\k\>p 


(2.13) 


Last  integrals  in  (2.13)  can  be  calculated  by  the  use  of  identity  (2.1 1).  In  (2.13)  Q,  P  and  K  are  the  values 
of  the  “losses”  of  the  constants  of  motion.  One  can  present  Snl  in  the  following  form 


Sni  (k)  =  Fk-  yknk 

(2.14) 

Fk  =  j1  Tkkfak,  \2S{cok  +  coki  +  coki  +  0)ki  )S(k+kl  +k2  +k3)x 
nknknkdkldk2dki 

(2.15) 

n  =  /l  Tkkikih  \2S{cok  +  coki  +  0)ki  +  o)kj  )S(k  +kl+k2+  k3) x 

K,n*2  +nknk3  ~nknh  )dkxdk2dk3 

(2.16) 

By  definition 


SHt)  =  F?,-ri,}n:  (2.17) 

F^p)  =  J I  Tkkik ^  I  8 {(Of.  +  Q)k  +  (Ok  +  (Ok  )S(k  +  iif  +  k2  +  k^x 

\k,\<p,\k2\<p,\k2\<p  (2.18) 

nknknkdkxdk2dk3 

rip)  =  j  I  T^hk 3  \28{(0k  +  G)ki  +  0)ki  +  0)k}  )8(k  +k{+k2+  k3)  x 

\ki\<p,\k2\<pt\k^\<p  (2.19) 

+\%  ~nknh)dkfdk2dk3 

One  should  mention  that  (k)  *  0  at  I  k  l>  p .  As  far  as  nk  =  0  at  I  k  l>  p ,  one  has: 

S„i(k)  =  Fkp,\k\>  p  (2.20) 

Formula  (2.20)  is  extremely  important.  It  expresses  the  following  clear  physical  fact:  the  income  term  Fk 
is  nonzero  far  beyond  the  domain  of  \k\<  p  where  the  wave  spectrum  is  concentrated.  Meanwhile 
Fk  ^  0  only  in  a  finite  domain.  Indeed,  vector  k  satisfies  the  conditions 


k  =  k2  +  k3-k} 

(2.22) 

(Ok  =cok2  +0)k3-0)ki 

(2.23) 

and  1  fc,  l<  p ,  1  k2  k  p  ,  1  k3  k  p . 

Conditions  (2.22),  (2.23)  can  be  satisfied  only  for 

1  *  k  Pmax 

(2.24) 

Here  pimx  =  f(p)  -  some  upper  limit  depending  on  a  shape  of  COk  . 
pmax  .  From  (2.22)  one  can  get 

.  One  can  get  some  apriori  estimate  for 

Pmax  <3/> 

(2.25) 

From  (2.23)  one  obtains 

co.  <  2  con 

P  max  P 

(2.26) 

More  accurate  estimate  for  pmM  is  defined  by  0)k  and  depends  on  the  depth.  On  the  infinite  depth 
COk  =  8^  a°d  P,mx  is  achieved  if 


(2.27) 


k2=k3,  kl  =  —  “ k2 1  \k2\=  p 


In  this  case 


(2.28) 

(2.29) 


Introducing  polar  coordinates  in  k  -pane,  we  have  the  following  expressions  for  losses  (rates  of  leakage)  of 
the  constants  of  motion  from  the  domain  I  k  fc>  p 


Pmax 

Q(P)=  J  pdpjFp(p,6)d0 

P  o 

Pmax  2jr 

P(P)=  J  Ptopdp  j Fp (p,0)d0 

P  o 

Pm&x  2tT 

Kx(p)  =  J  p2dp j Fp(p,0)cos0d0 

P  0 

Pmax  2/r 

Ky(p)=  J  p2dpjFp(p,0)sin0d0 

p  o 


(2.30) 


3.  Clean  test  for  integrals  conservation 


Now  we  can  answer  the  question  about  real  conservation  of  the  integrals.  If  the  domain  is  finite,  they  are 
never  preserved.  It  is  obvious  from  (2.30)  that  in  all  cases 

Q(p)>0,  P(p)> 0  (3.i) 

Thus,  wave  action  and  energy  always  leak  out  of  the  domain  I  k  l<  p  .  In  most  cases,  wave  spectrum  is 

K  71 

concentrated  in  the  right  half-plane  -  —  <  0  <  — ,  COS 0  >  0 .  In  this  situation  Kx{p)>  0 ,  the  sign  of 

Ky(p)  can  be  arbitrary. 

Suppose  again  that 

nk  =0,\k\>  p  (3.2) 


We  showed  that 

^  =  S(np\k)  =  F{p)(k)>  0 


(3.3) 


if  \k\<Pma. 


Hence,  from  physical  viewpoint  condition  (3.2)  is  artificial.  If  it  is  satisfied  in  the  initial  moment  of  time,  it 
will  be  immediately  violated.  In  the  close  moment  of  time  dt 


n  =  n0  +  S£,(k)  ■  dt 


(3.4) 


n  becomes  positive  in  the  whole  domain  I  k  k  pmax . 

Anyway,  the  consideration  we  performed  is  useful.  It  can  be  used  as  a  foundation  for  a  “clean”  test 
for  all  codes  for  numerical  solution  of  the  kinetic  wave  equation. 

If  condition  (3.2)  is  satisfied 

dn  _  , ,  . 

—  =  0  ii\k\>pmm  (3.5) 


Hence 


a_ 

dt 


l*l<Pm 


jsidk 

'P  max 

jt  j<»tsZdk 


=  0 


=  0 


(3.6) 


r=0 


dt 


\ks;,dk 


=  0 


(=0 


Conditions  (3.6)  can  be  rewritten  in  polar  coordinates  as 

Pmax  2tT 


-  J  pdpjs'fde 


0 

Pmax  'Ik 

-  j  Q)ppdp\s{fd0 


=  0 


\t= 0 


dt 

d_ 

dt 

d_ 

dt 


=  0 


l/=0 


0 

Pmax  2  K 

J  p2dp  j S(J}  cosddd 


0 

0  Pinax 
0 


0 

2  K 


(3.7) 


=  0 


1 1=0 


—  j  p2 dp  j  S{f  sin  0d6 


=  0 


1=0 


Condition  (3.7)  can  be  relatively  easily  checked  for  any  numerical  code  used  for  solution  of  the  kinetic 
equation.  To  check  the  quality  of  the  code  one  should  put  the  initial  data  n.  >  0  at  I  A'  k  n 
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expansion  of  the  integration  domain  is  the  price  to  be  paid  for  nonlocality  of  the  four-wave  interactions. 


Relations  (2.13)  can  be  generalized  to  the  case  when  there  is  the  interaction  with  wind  and 
damping.  Now  kinetic  equation  (2.1)-(2.2)  reads 


^~Snt+  &n* 


(3.9) 


where  is  growth-rate  of  the  instability  or  the  damping  depending  on  the  sign.  Equations  (2.6)-(2.8) 
now  should  be  replaced  by  the  following  relations: 


(3.10) 


Relations  (3.10)  are  formal  and  for  a  finite  domain  should  be  deciphered  in  a  proper  way.  To  do  that  one 
can  assume 


Pk  =  -A ,  A  — »  ,  |  k  l>  p 

In  the  domain  I  k  \>  p  one  can  neglect  the  time  derivative  —  and  put 

dt 

S„-Plnl=F'-{yt+A)nt^  0 

As  far  as  yk  «  A ,  one  can  consider  approximately  ( I  k  l>  1 ) 


By  plugging  (3.1 1)  into  (3.10)  one  notice  that  the  value  of  A  is  cancelled  from  the  equations,  taking  the 
form 


dNp 

dt 

dEp 

dt 

dM 

dt 


-  j  Pknkdk  -  Q 

\k\<p 

~  \<°kPknkdk-P 

\k\<p 

-  jkfiknkdk  -  K 

\k\<p 


(3.12) 


Equations  (3.12)  are  the  balance  equations  which  also  could  be  used  to  control  numerical  codes.  One  should 
either  specially  program  the  calculation  of  the  losses  Q,P,K  or  extend  the  integration  domain  to 

I  k  1=  • 

i  max 

4.  Estimates  for  integral  losses. 


Suppose  that  spectral  density  of  wave  action  n{k)  has  the  maximum  at  k  ~  k0  and  p»  k0.  Let  us 

estimate  in  this  limit  the  losses  Q,  P  and  K.  Integrals  (2.30)  consist  of  two  parts.  One  part  is  given  by 
integration  in  a  small  domain  p'=p  +  8p,  dp  <  k0 .  In  this  domain,  integration  by  k{ ,  k2  is  performed 
over  the  vicinity  of  the  spectral  maximum  k0 .  Thus 


e=e(1)+e(2) 

p _ pm  +pM 

—  —  0)  —  (2) 

K  =  K  +K 


Here 


Qm  r  ”Q)»  (*o) 
co(k0) 


PtO’Pto 


Pm  *^n(P)n\ka)\T  l!'.‘ 


Q)(k0) 

P 

o)(k0) 


o '  I  P^Ph 


Xa)  = — — — n(p)n2(k)\T 


P<k0,p,k0 


k6 

K o 


(4.1) 


(4.2) 


Other  part  of  contribution  in  (4.1)  is  given  by  integration  in  the  domain  p'  =  p  An  this  case,  all  wave 
vectors  in  (2.18)  have  the  same  order  of  magnitude.  Hence 


i(2)  __  n  (p) 


Qw  =■ 

0>(p) 

Pm=n\p)\Trrrpfp‘ 


P>P>P,P\ 
|3 


co(p) 


n\p)-T  2  p6 

P,P*P*P  V 


We  consider  now  the  case  of  deep  fluid.  In  this  case  0)(p)  =  g^2 p^2 
T(P’P>P,P)  =  P3 

„  kn  «  p 
T(p,k0,p,k0)  =  k^p 


(4.3) 


(4.4) 


We  will  consider  that  n(k )  is  a  powerlike  function 


n(k)  =  k  5 


(4.5) 


Comparing  Q(V>  and  Qa)  one  see  that  for  the  case  S  >  19/4  Qm  «  Q(2> .  The  same  is  correct  for 
other  constants  of  motion.  Finally,  one  obtains 

Q  =  p23'2^ 

P  =  pl2~3s  (4.6) 

k  =  p25/2-3s 

Now  we  can  answer  the  question  about  real  conservation  of  the  motion  constants.  All  motion  constants 
conserve  if 

n{k)  <  Ck~25^6  at  k  — >  °°  (4.7) 


5.  Kolmogorov  spectra  and  their  experimental  confirmation. 


It  is  known  since  1966  (Zakharov,  Filonenko)  that  the  stationary  equation 

S„=0 

has  isotropic  powerlike  solutions 

n?=C,k-2  3'6 

n(2)  =  C2k-* 


(5.1) 


(5.2) 

(5.3) 


The  physical  meaning  of  this  solution  becomes  clear  after  plugging  (5.2)  into  (4.6).  For  s  =  23/6  Q  is 
the  constant  while  P  and  K grow  in  time.  Hence,  (5.2)  is  a  Kolmogorov  spectrum,  corresponding  to  constant 

flux  of  wave  action  from  large  to  small  wavenumbers.  For  such  type  of  asymptotic  neither  N ,  E  or  M  are 
“real”  constants  of  motion. 

If  J  =  4  Q  =  p  1//2  .  In  this  case  Q  — >  0  at  p  — >  °o  ,  and  wave  action  is  “a  real”  constant  of 
motion.  Meanwhile,  in  this  case  P  —  const .  Hence,  (5.3)  is  a  Kolmogorov  spectrum  describing  permanent 
leakage  of  energy  to  large  wavenumbers.  In  terms  of  spectral  density  of  energy  spectra  (5.2)-(5.3)  read 

C^e^V"'3  (5.4) 

E™  =  a2Pi/3a)-4  (5  5) 

Here  al ,  a2  are  unknown  Kolmogorov  constants.  One  can  formulate  a  conjecture  that  a  general 
physically  relevant  solution  equation 

nk  =  p*5t-‘F(£(*i)1«,  j(gk)-"2) 


(5.6) 


where  F  is  an  unknown  function  of  two  variables.  It  can  be  found  explicitly  in  heuristic  “diffusive”  model  of 
Snl  (Zakharov,  Pushkarev  1999). 
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The  spectrum  CD  has  the  long  history.  Zakharov  and  Filonenko  found  it  analytically  as  a  well- 
hidden  exact  solution  of  equation  (5.1)  in  1966.  Both  authors  lived  in  the  USSR  and  were  not  allowed  to 
travel  abroad  and  report  their  results  on  international  conferences.  That  was  one  reason  why  the  paper  of 
Zakharov  and  Filonenko,  published  in  the  leading  Russian  scientific  journal  (Doklady  Academii  Nauk)  was 
almost  not  noticed. 

There  was  another  reason  too.  In  1958,  O.Phillips  offered  that  the  spectrum  of  wind-driven  waves 
is  defined  completely  by  wave  breaking  and  has  a  universal  form 

e.-ogV5  (5.7) 

Here  (X  is  the  dimensionless  “Phillips  constant”. 

The  very  idea  of  Phillips  was  seminal  and  productive.  A  we  understand  now,  situations  when 
spectra  of  wave  turbulence  are  defined  completely  by  local  singularities  is  rather  common  (see,  for  instance, 
Zakharov,  Dias,  Pushkarev,  Guyenne  2000). 

According  to  original  idea  of  Phillips,  spectrum  (5.7)  is  automatically  established  at  CO>  O)0, 
where  Q)0  is  a  characteristic  frequency  of  the  spectral  maximum.  All  spectral  dynamics  is  just  evolution 

=  t»0(0  0  at  f  .  ft,  spite  of  its  elegance  and  simplicity,  the  initial  conjecture  of  Phillips  is  not 

confirmed  by  experiments.  By  definition 

oo 

\eJo)  =  {n2)  =  H2  (5.8) 

0 

By  integration  of  (5.7)  one  obtain 

H2  ~  a$2  s  aS2  r4 
4  6) 4  4(2;r)4 


or 


Here  T  —  2 7cj COq  -  characteristic  period  of  energy  containing  waves. 

Relation  (5.9)  can  be  checked  experimentally.  It  was  Yoshiaki  Toba  who  did  it.  In  1972-73,  he 
published  a  set  of  articles,  summarizing  his  long-time  experiments  on  the  wind-wave  channel  in  Sendai 
University.  He  found  that,  instead  of  (5.9),  another  relation  holds: 


HocT3/2 


(5.10) 


Moreover,  the  Phillips  constant,  supposed  to  be  universal  and  not  depending  on  the  wind  velocity,  happened 
to  be  proportional  to  the  “  friction  velocity”  ut ,  characterizing  the  momentum  transfer  from  air  to  sea. 
These  facts  can  be  naturally  explained  if  one  assume  that  instead  of  (5.7),  the  spectrum  has  a  form 


eco  = 


(5.11) 


Here  fi  is  another  dimensionless  constant  (Toba’s  constant).  Toba  made  careful  measurements  of  the 
spectrum  tale  and  found  that  formula  (5.9)  describes  the  spectral  asymptotic  very  well.  He  found 
experimentally  that  the  value  of  ($  is: 

/?  =  6.2  •  10-2  (5.12) 

We  must  stress  that  Toba  was  completely  unaware  about  the  paper  of  Zakharov  and  Filonenko. 

After  pioneering  works  of  Toba  experimental  confirmation  of  CO spectrum  are  mounted.  Just  a 
list  of  authors  who  observed  this  asymptotic  is  impressive:  among  them  Kawai  at  al.  1977,  Kahma  1981, 
Forristall  1981,  Battjes  at  al  1986,  Donelan,  Hamilton  and  Hui  1985,  Donelan  and  Pierson  1987.  In  1985 
O.Phillips,  summarizing  all  experimental  data,  criticized  his  early  theory  (which,  in  our  opinion,  is  still  a 

sample  of  outstanding  scientific  intuition)  and  agreed  that  G !/  is  a  reality. 

Another  approach  for  confirmation  of  the  Of  4  spectrum  is  the  numerical  experiment.  Several 
groups  (see  Masuda  and  Komatsu  1980,  1996;  Resio  and  Perrie  1991;  Polnikov  1994,  2000)  observed  a 

universal  effect.  If  the  initial  data  in  equation  (2.1)-(2.2)  decay  like  G)~5 ,  very  soon  its  asymptotic  behavior 
changes  to  Of4 . 

In  this  connection,  we  would  like  to  mention  especially  works  of  the  group  of  D. Resio.  They  did 
not  only  follow  formation  of  the  Ctf*  spectrum,  but  also  calculated  fluxes  of  energy  at  k  — >  °°  and 
checked  that  for  Of4  the  flux  is  constant  in  G) . 

In  conclusion,  one  can  say  that  the  spectrum  G)  4  is  definitely  confirmed  now  by  both  many 

experiments  as  well  as  many  numerical  simulation  arguments.  The  spectrum  0)~nn  also  is  confirmed  quite 
well,  but  this  point  is  beyond  the  scope  of  this  article. 


5.  Why  of4  not  ct)~5  ? 


The  battle  between  Of4  and  G)  5  could  look  strange  for  a  person  outside  of  a  narrow  community  of 
Snl  experts.  Nevertheless,  this  argument  makes  a  serious  sense.  It  is  enough  to  compare  the  value  of  the 

energy  loss  P  in  (4.6)  for  both  spectra.  For  £a  =  Of  4  and  s  =  4 ,  P  is  the  constant.  .  For  £a  =  Of  5 
5  =  9/4,  P  =  p~3/2  — >  0  at  p^>  °o . 

Another  words,  for  £a  =G)  5  the  energy  in  the  wave  ensemble  is  conserved,  while  in  the  case 

—4 

£a  =  G)  it  leaks  out  with  a  constant  rate.  This  is  a  very  critical  difference.  Just  an  elementary  analysis  of 

the  observational  data  contradicts  the  idea  that  both  wave  action  and  energy  are  conserved.  Indeed,  all 
experiments  show  that  the  spectral  maximum  moves  in  the  process  of  “maturing  “  to  the  low  wave  numbers. 
At  the  same  time,  the  spectrum  in  the  asymptotic  area  G)  »  G) q  stays  almost  constant.  As  far  as  quanta  of 

waves  loose  their  energy,  moving  from  high  to  low  frequency  region,  the  outlined  facts  are  compatible  with 
the  fact  of  permanent  loss  of  wave  energy,  existence  of  the  constant  flux  of  energy  to  high  G)  and,  as  a 
result,  to  the  Of4  asymptotic  at  G)  — >  oo  . 

Summarizing  the  facts,  one  can  say  that  the  asymptotic  £a  =  Of5  in  the  weak  turbulent  regime  is 
a  contradiction  to  the  energy  conservation  law. 

Said  Aristotle:  “You  are  friend  Plato,  but  the  truth  is  more  valuable”.  We  can  say 
You  are  a  friend  Dr....  but  the  conservation  of  energy  is  more  important”. 
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Abstract 

A  two-parameter  nonlinear  dispersive  wave  equation  proposed  by  Majda,  McLaughlin 
and  Tabak  is  studied  analytically  and  numerically  as  a  model  for  the  study  of  wave  turbu¬ 
lence  in  one-dimensional  systems.  Our  ultimate  goal  is  to  test  the  validity  of  weak  turbulence 
theory.  Although  weak  turbulence  theory  is  independent  on  the  sign  of  the  nonlinearity  of 
the  model,  the  numerical  results  show  a  strong  dependence  on  the  sign  of  the  nonlinearity. 

A  possible  explanation  for  this  discrepancy  is  the  strong  influence  of  coherent  structures  - 
wave  collapses  and  quasisolitons  -  in  wave  turbulence. 

1  Introduction 

A  wide  variety  of  physical  problems  involve  random  nonlinear  dispersive  waves.  The  most 
common  tool  for  the  statistical  description  of  these  waves  is  a  kinetic  equation  for  squared 
wave  amplitudes,  the  so-called  kinetic  wave  equation.  Sometimes  this  equation  is  also  called 
Boltzmann’s  equation.  This  terminology  is  in  fact  misleading  because  the  kinetic  wave  equation 
and  Boltzmann’s  equation  are  the  opposite  limiting  cases  of  a  more  general  kinetic  equation 
for  particles  which  obey  Bose-Einstein  statistics  like  photons  in  stellar  atmospheres  or  phonons 
in  liquid  helium.  It  was  first  derived  by  Peierls  in  1929  [1].  In  spite  of  the  fact  that  both 
the  kinetic  wave  equation  and  Boltzmann’s  equation  can  be  derived  from  the  quantum  kinetic 
equation,  the  kinetic  wave  equation  was  derived  independently  and  almost  simultaneously  in 
plasma  physics  and  for  surface  waves  on  deep  water.  This  was  done  in  the  early  sixties  while 
Boltzmann’s  equation  was  derived  in  the  nineteenth  century!  The  derivation  for  surface  waves 
is  due  to  Hasselmann  [2]  (see  also  Zakharov  [3]). 

Once  the  kinetic  wave  equation  has  been  derived,  the  shape  of  wave  number  spectra  can 
be  predicted  by  the  so-called  weak  turbulence  (WT)  theory.  It  is  called  weak  because  it  deals 
with  resonant  interactions  between  small-amplitude  waves.  Thus,  contrary  to  fully  developed 
turbulence,  it  leads  to  explicit  analytical  solutions  provided  some  assumptions  are  made.  So 
far,  there  have  been  only  a  few  studies  to  check  the  results  of  WT  theory.  Recently,  Pushkarev 
and  Zakharov  [4]  numerically  solved  the  three-dimensional  dynamical  equations  for  the  free- 
surface  elevation  and  the  velocity  potential  in  the  case  of  capillary  water  waves.  They  obtained 
an  isotropic  spectrum  close  to  the  theoretical  power-law  found  by  Zakharov  and  Filonenko 
[5].  Majda,  McLaughlin  and  Tabak  [6]  (hereafter  referred  to  as  MMT)  considered  four- wave 
interactions  by  introducing  a  one-dimensional  model  equation.  This  equation  can  be  integrated 
numerically  quite  efficiently  on  large  inertial  intervals.  They  examined  a  family  of  Kolmogorov- 
type  solutions  depending  on  the  parameters  of  the  equation.  The  validity  of  several  theoretical 
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hypotheses  was  then  assessed  numerically.  Namely,  MMT  confirmed  the  random  phase  and 
quasi-gaussian  approximations.  They  also  showed  the  independence  of  the  solutions  on  the 
nature  of  forces,  initial  conditions,  and  the  size  and  level  of  discreteness  of  the  computational 
domain.  However,  their  simulations  surprisingly  displayed  spectra  steeper  than  the  predicted 
ones.  They  explained  the  discrepancy  by  proposing  a  new  inertial  range  scaling  technique  which 
seems  to  yield  the  appropriate  exponents.  More  recently,  Cai,  Majda,  McLaughlin  and  Tabak 
[7]  revisited  their  earlier  results  and  found  some  results  which  agree  with  WT  theory  as  well.1 
They  considered  two  kinds  of  Hamiltonians:  Hamiltonians  which  are  the  sum  of  a  quadratic 
term  and  a  quartic  term  (positive  nonlinearity),  as  in  [6],  and  Hamiltonians  which  are  the 
difference  between  a  quadratic  and  a  quartic  term  (negative  nonlinearity).  For  either  sign  of 
nonlinearity,  they  found  agreement  with  MMT  theory  in  some  cases  and  agreement  with  WT 
theory  in  some  other  cases.  Since  their  computations  were  performed  with  a  dispersion  relation 
in  which  the  frequency  varies  like  the  square  root  of  the  wave  number,  one  can  see  an  analogy 
with  deep  water  waves.  Incidentally,  the  WT  theory  was  recently  developed  for  shallow  water 
waves  by  Zakharov  [8]. 

As  in  many  other  fields,  numerical  modeling  leads  to  some  difficulties,  especially  when  one 
wants  to  compare  with  the  theory.  Most  of  these  difficulties  are  related  to  finite-size  effects,  i.e. 
the  domain  is  discretized  into  a  grid  of  points  in  computations  whereas  one  assumes  an  infinite 
medium  in  theory.  We  can  mention  the  bottleneck  phenomenon  [9]  which  tends  to  flatten  the 
slope  of  the  inertial  range  at  small  scales.  It  is  commonly  observed  in  problems  with  a  dissipative 
cutoff.  In  addition,  Pushkarev  [10]  revealed  the  phenomenon  of  frozen  turbulence  at  very  low 
levels  of  nonlinearity.  In  this  situation,  the  resonance  conditions  have  very  few  solutions  (or 
may  not  be  fulfilled  at  all!)  because  of  the  discrete  values  of  wave  numbers.  As  a  consequence, 
there  is  no  energy  flux  due  to  the  lack  of  resonating  wave  vectors.  The  power-law  regime  only 
takes  place  at  moderate  levels  of  nonlinearity  where  quasi-resonant  interactions  come  into  play. 
Pushkarev  concluded  that  weak  turbulence  in  bounded  systems  combines  the  features  of  both 
frozen  and  Kolmogorov-type  turbulence.  The  beauty  of  the  MMT  model  equation  is  that  the 
above  mentioned  difficulties  can  be  controlled  completely. 

After  introducing  the  model  equation,  the  paper  is  divided  into  two  parts.  In  the  first 
part,  the  MMT  equation  is  studied  analytically.  A  weak  turbulence  description  of  the  equation 
is  provided  (see  [6]).  We  find  the  Kolmogorov  solutions  of  the  kinetic  equation  and  determine 
the  set  of  parameters  for  which  such  solutions  can  be  realized.  Then  we  discuss  the  coherent 
structures  which  can  compete  with  weak  turbulence.  The  most  simple  coherent  structures  are 
solitons  similar  to  the  soliton  solutions  of  the  Nonlinear  Schrodinger  Equation  (NLS). 

Solitons  for  the  MMT  equation  exist  only  if  nonlinearity  is  negative.  In  the  cases  of 
interest,  they  are  shown  to  be  unstable  (see  Section  7)  and  cannot  play  an  important  role  in 
the  wave  dynamics. 

As  an  alternative  to  soliton  coherent  structures,  there  are  wave  collapses  described  by 
self-similar  solutions  of  the  MMT  equation.  These  solutions  can  exist  in  a  certain  parameter 
regime  for  both  signs  of  nonlinearity.  Theoretically  speaking,  both  solitons  and  collapses  can 
coexist  with  weak  turbulence. 

Another  type  of  coherent  structures  are  quasisolitons,  or  envelope  solitons.  They  were 
discussed  recently  by  Zakharov  and  Kuznetsov  [11].  In  the  MMT  model  quasisolitons  exist  at 
positive  nonlinearity  only.  Their  stability  remains  an  open  question. 

The  main  new  theoretical  results  of  the  first  part  are  a  careful  tabulation  of  the  signs 

1  These  three  papers  were  kindly  given  to  us  when  the  present  manuscript  was  essentially  completed.  Some 
of  the  results  are  similar  to  ours,  but  their  interpretation  is  different. 
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of  the  fluxes  for  the  MMT  model  equation,  the  existence  and  possible  role  of  quasisolitons  for 
positive  nonlinearity,  and  an  analogy  with  Phillips  spectrum  associated  with  the  formation  of 
collapses. 

In  the  second  part,  we  describe  the  results  of  the  numerical  study  of  the  MMT  equation. 
We  find  that  the  wave  turbulence  described  by  the  MMT  equation  is  different  both  quanti¬ 
tatively  and  qualitatively  for  both  signs  of  nonlinearity.  Since  the  predictions  of  WT  theory 
are  identical  for  both  signs  of  nonlinearity,  WT  theory  can  be  applied  at  best  for  one  sign  of 
nonlinearity.  Our  analysis  of  the  results  leads  to  somewhat  contradictory  results. 

For  positive  nonlinearity  the  balance  of  energy  and  particle  fluxes  as  well  as  the  level  of 
turbulence  are  in  good  agreement  with  WT  theory.  Meanwhile  the  slope  of  the  spectrum  in  the 
window  of  transparency  is  steeper  than  predicted  by  WT  theory. 

In  the  case  of  negative  nonlinearity  the  picture  of  turbulence  is  quite  different  from  the 
WT  predictions,  both  qualitatively  and  quantitatively.  First  of  all,  the  turbulence  is  stabilized 
on  a  level  which  is  one  order  of  magnitude  less  than  predicted  by  WT  theory.  Then  the  sign 
of  the  flux  of  particles  is  opposite  to  the  one  predicted  by  WT  theory.  Both  these  facts  lead  to 
a  conjecture  on  the  existence  of  a  strong  and  essentially  nonlinear  mechanism  which  competes 
successfully  with  WT  quartic  resonances.  In  our  opinion,  this  mechanism  is  the  wave  collapse, 
described  by  self-similar  solutions  of  the  MMT  equations.  At  the  same  time,  the  high-frequency 
tail  of  the  spectrum  has  a  slope  which  coincides  exactly  with  the  slope  predicted  by  WT  theory. 
This  leads  to  the  conclusion  that  in  spite  of  the  presence  of  wave  collapses,  the  high-frequency 
asymptotics  of  spectra  is  governed  by  the  WT  processes  which  are  responsible  for  carrying  only 
a  small  part  of  the  energy.  The  coexistence  of  wave  collapses  and  weak  turbulence  was  already 
described  in  the  context  of  the  2D  NLS  [12]. 

Wave  collapse  is  an  example  of  an  essentially  nonlinear  coherent  structure  arising  in 
wave  turbulence  under  certain  conditions.  As  said  above,  another  important  type  of  coherent 
structures  are  quasisolitons  or  envelope  solitons  living  for  a  finite  time.  Such  structures  can 
arise  in  the  MMT  model  in  the  case  of  positive  nonlinearity.  We  believe  that  these  structures 
are  responsible  for  the  deviation  of  the  spectra  from  the  ones  predicted  by  WT  theory. 


2  Model  equation 


We  investigate  the  family  of  dynamical  equations 
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where  ^{x.t)  denotes  a  complex  wave  field  and  a,  /3  are  real  parameters. 

If  A  —  +1,  one  exactly  recovers  the  MMT  model  which  was  treated  in  [6].  Note  that 
our  parameter  (3  is  the  opposite  of  the  parameter  /?  in  MMT.  The  extension  A  =  ±1  in  Eq. 
(2.1),  which  was  also  treated  in  [7],  raises  an  interesting  problem  because  the  balance  between 
nonlinear  and  dispersive  effects  may  change  according  to  A. 

Besides  the  Hamiltonian 


H  —  Hl  +  Hnl  = 


/ 


d_ 

dx 


dx , 


(2.2) 
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the  system  (2.1)  preserves  two  other  integrals  of  motion:  wave  action  and  momentum,  respec¬ 
tively 

N=j\*U*  and  ^=j/(+§£  -Yxr)ix- 

As  usual,  it  is  convenient  to  work  in  Fourier  space.  Let  us  write  Eq.  (2.1)  as 

i  =  w(k)i>k  +  J  Tnzk  i’lfoi’z  +  k2~k3~k)  dkidk2dk3 ,  (2.3) 

where  ipk  =  V>(M)  denotes  the  fc-th  component  in  the  Fourier  decomposition  of  i>(x,t)  and  (*) 
stands  for  complex  conjugation. 

In  this  form,  Eq.  (2.3)  looks  like  the  so-called  one-dimensional  Zakharov’s  equation 
determined  by  the  linear  dispersion  relation 

cj(k)  =  l/cl01  ,  a  >  0,  (2.4) 

and  the  simple  interaction  coefficient 

Tmfe  =  T(ku k2, k3,k)  =  A  Ihkihkf*  .  (2.5) 

One  easily  sees  that  the  kernel  Ti23fc  possesses  the  symmetry  required  by  the  Hamiltonian 
property 

?123 k  =  ^213 k  =  ^12fc3  =  ?3fcl2  •  (2-6) 

Moreover,  the  absolute  values  in  Eqs  (2.4)  and  (2.5)  ensure  the  basic  assumptions  of  isotropy 
and  scale  invariance.  In  other  words,  a j(k)  and  Ti23fc  are  invariant  with  respect  to  rotations 
(k  — >  —k)  and  they  are  homogeneous  functions  of  their  arguments  with  degrees  a  and  /3 
respectively,  i.e. 

T(tkutk2,tk3,Zk)=Zl3T(kuk2,k3,k),  <£>0.  (2.7) 

Following  MMT,  we  fix  a  =  1/2  by  analogy  with  gravity  waves  whose  dispersion  relation  reads 
as  uj  =  (gk)1/2  ( g  being  the  acceleration  due  to  gravity).  The  power  (3  takes  the  value  +3  if 
the  analogy  is  extended  to  the  nonlinear  term  but  we  will  consider  a  wider  range  of  values  for 
/3. 

Eq.  (2.3)  describes  four-wave  interaction  processes  obeying  the  resonant  conditions 

k\  4-  fc2  =  &3  +  k ,  (2-8) 

4“  0J2  =  Ct?3  + 1 o  .  (2.9) 

For  a  >  1  these  equations  only  have  the  trivial  solution  £3  =  fci,  k  =  fc2  or  £3  =  fc2,  k  =  k\.  For 
a  <  1  there  is  also  a  nontrivial  solution.  Note  that  in  this  case  the  signs  of  ki  must  be  different. 
For  instance,  k\  <  0  and  fc2,  ks^k  >0.  If  a  =  1/2,  Eqs  (2.8)-(2.9)  can  be  parametrized  by  two 
parameters  A  and  £ 

h  =  -A2e,  k2  =  A2(  l  +  £  +  £2)2,  h  =  A2{  l+£)2,  k  —  A2£2(l  +  £)2 .  (2.10) 

In  the  case  a  =  2  and  /?  =  0,  Eq.  (2.1)  becomes  the  NLS  equation 

*£-S  +  *W*  P.n, 

(note  here  that  |^|2  =  -^r)- 

Positive  nonlinearity  A  =  +1  corresponds  to  the  defocusing  NLS,  while  negative  nonlin¬ 
earity  corresponds  to  the  focusing  NLS. 
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PART  I:  THEORY 


i 


3  Weak  turbulence  description  of  the  model  equation 

If  one  only  considers  small  nonlinear  effects,  then  the  statistical  behavior  can  be  mainly  de¬ 
scribed  by  the  evolution  of  the  two-point  correlation  function 

(^(fc,  t)^*(A;',t))  =  n(k ,  t)  8(k  -  kf) , 

where  brackets  denote  ensemble  averaging.  We  introduce  also  the  four- wave  correlation  function 

${ki,t)i>(h2,t)i>*(k3,t)i>*{k,t))  =  JuzkKh  +k2-k3-k).  (3.1) 

On  this  basis,  WT  theory  leads  to  the  kinetic  equation  for  n(fc,  t)  and  provides  tools  for 
finding  stationary  power- law  solutions.  For  details,  see  [6].  Here  we  explain  the  main  steps  of 
the  procedure  applied  to  our  model. 

The  starting  point  is  the  original  equation  for  n(fc,t).  From  Eq.  (2.3),  we  have 

=  2  J  Im  Jyizk  Ti23k  tf(*i  +  k2  -  k3  -  k)  dkidk2dk3 .  (3.2) 

Due  to  the  Quasi-Gaussian  Random  Phase  approximation 

Re  Ji23fc  -  riin2  [8(ki  -  fc3)  +  8(k\  -  k)]  .  (3.3) 

The  imaginary  part  of  J\2%k  can  be  found  through  an  approximate  solution  of  the  equation 
imposed  on  this  correlator.  The  result  is  (see  [13]) 

Im  J±2Sk  —  27tT1*23 k  +  ^2  —  ^3  “  w)  (711112713  +  n\ri2nk  —  niu^rik  —  ft2n3nfc) .  (3.4) 

This  gives 

|Ti23*|2  (™in2™3  +  niU2Uk  -  nxnzUk  -  n2n3nfc)  (3.5) 

x  5(w i  +  U2  —  cos  —  oj)  S(k\  +  fc2  “  fc3  —  k)  dk\dk2dk^  . 

Since  the  square  norm  cancels  the  sign  of  T123a;,  it  is  clear  that  the  WT  approach  is  independent 
on  A.  Here  we  point  out  that  MMT  mistakenly  wrote  a  factor  127 r  instead  of  47 r  in  Eq.  (3.5) 
and  the  right  hand  side  of  Eq.  (3.5)  with  the  opposite  sign.  This  fact  is  particularly  important 
when  determining  the  fluxes  of  wave  action  and  energy. 

Assuming  that  n(-k)  =  n(k)  (similarly  to  an  angle  averaging  in  higher  dimensions),  one 

gets 


diik  A  [ 

Ht  =  4'/ 


W(»>) 

dt 

x 

+ 

+ 


(nin2n3  +  nin2nw  -  nin3nw  -  n27i3nw) 

+  CJ2  —  td3  —  Cd)  ^(tdi  <*  +  U)2  —  Cd3  o  +  id  a  ) 

r  /  —  .  X  J.  X\  i  XXX  _1_  , 

a  +  w2«  +  a;3  a  —  a; « )  +  o(a;i  <*  —a;2a  —  a;3  ) 

6(— wi«  +W2“  —  —  w«)|  dw\dw2du3, 


in  f 

^4  /  (W1W2W3W)  <* 


(3.6) 
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where  Af{ o>)  —  n  (k(cj))  dk/du ;,  stands  for  n  (k(uj))  and  u  is  given  by  Eq.  (2.4). 

The  next  step  consists  in  inserting  the  power-law  ansatz 

n(  a/)  oca;-7,  (3.7) 

and  then  performing  the  Zakharov’s  conformal  transformations  [6,  12,  13].  Finally,  the  kinetic 
equation  becomes 

ccu^-1 7),  (3.8) 

where 

Ja 

x  5(£i«  +  £2“  +  £3“  “  1)  (1  +  £3^  -  £iy  -  £2^)  d£id£2d£3  ,  (3.9) 


with 


A  =  {0  <  fi  <  1 , 0  <  &  <  1 ,  £1  +  £2  >  1}  and  y  =  37  +  1  -  —  +  3  . 

a 

The  nondimensionalized  integral  /(a, /3,7)  is  obtained  by  using  the  change  of  variables  c Oj  — ► 

(j  =  1, 2, 3) . 

The  ansatz  (3.7)  makes  sense  if  the  integral  in  (3.6)  converges.  It  could  diverge  both  at 
low  and  high  frequencies.  The  condition  of  convergence  at  low  frequencies  coincides  with  the 
condition  of  convergence  of  the  integral  in  (3.9)  and  can  be  easily  found.  It  reads 


27  <  -1  + 


/?  +  4 
a 


(3.10) 


The  condition  of  convergence  at  high  frequencies  can  be  found  after  substituting  (3.7)  into  (3.6). 
Omitting  the  details,  we  get  the  result 


^  P  +  a- 1 

7  >  - 

a 


(3.11) 


In  all  the  cases  discussed  in  this  article,  both  conditions  (3.10)  and  (3.11)  are  satisfied. 
For  the  case  a  —  1/2,  one  can  transform  Eq.  (3.6)  into  the  form 


=  647r  w4(/9+l)  [5l  +S2  +  S3  +  54]  , 
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where 


Si 


-  »/ 


S2 


S3 


S4  = 


^^(u  +  l)2^1 


(r 


“  (rr^H?") n  (1 


u(u  +  1) 
+  u  +  u2 


+  n(cj)n 
-  n  (w)  n 


(1  +  «  +  w2)3^+4 

ti(u  +  1)  \  /  U  +  1  \  (  \  ( 

u(u  +  1) 


a;  n 


u 


+  u  +  u2 
u 

+  u  +  u2 


+  u  +  u2 


(  u  +  1 

\  1  +  U  +  U2  W 
\  /  ti+1 

win  - - 

J  \1  +  u  + 


7W 


-l 

+  n  (1 


1  (l-hti4-uY+1(u  +  l)2W 


u 


2/0+3 


/^1+u  +  u2  ^  „  ,  1X  \  /lfc  +  1  \ 

n  - (j  n  ((u  -f-  l)a;)  n  - w 

.  V  u  J  \  u  J 


//  1\  \  (u  +  1  \  ,  ,  / 1 -f  u  +  u2  \  /  w  +  1  \ 

w)n((u  +  l)w)n  ( — — cjJ  —  n(a/)n  I  - - - w  J  n  I  — — w  J 

-  n  (a;)  n  ^  ^ ^  n  ((u  +  l)w)  du , 

71  w2^+2(l  +u  +  w2)/3+1  f  ,  .  /l  +  u  +  u2  \  /  u  \ 

y0  (1  +  ^p+3  L  V  !  +  w  /  \i  +  «  / 

,  x  fl  +  u  +  u2  \  (  U  \  ( l  +U  +  V? 

(a?)  n  — - - u)  n  - - w  -  n  (u/J  n  f  uuj J  n  — - - 

v  '  \  1+u  J  \l+u  J  v  }  K  }  V  l  +  « 

n  (w)  n  (uw)  n  ^  j  du , 

fl  (1  +  u  +  u2)^+1  T  /w\  /1  +  u  +  u2  \  / 

y0  u2^+3(i  +  w)2^+2  [  \u/  \  u(i  +  w)  y  n  Vi 


+  n 


■cj 


UJ 


+  u 


u2^+3(l  + 

,  ,  \  fl+U  +  U2  \  fUJ\  {U\  (  w  V 

+  n  (an  n  — 7- - r-w  n  —  -  n  (w)  n  [  —  n  - 

v  }  \  u(l  +  u)  J  \uJ  v  }  \u)  Vl+uy 


-w 


—  n(u)n 


&) n  (L 


+  u  +  u4 


n  (w)  n  | 

-u>^j  du . 


u(l  +  u) 

This  equation  can  be  used  for  the  numerical  simulation  of  weak  turbulence. 


4  Kolmogorov  solutions 

The  aim  is  to  look  for  stationary  solutions  of  the  kinetic  equation.  From  Eq.  (3.8)  we  easily 
find  that  the  solution 


dJ\f  (w) 

Qt  =  0  I{a,fi,i)=  0,  (4.1) 

is  satisfied  only  for  7  =  0, 1  and  y  =  0,1. 

In  terms  borrowed  from  statistical  mechanics,  the  cases  7  =  0, 1  represent  the  thermody¬ 
namic  equilibrium  solutions 


n(w)  =  c, 

(4.2) 

where  c  is  an  arbitrary  constant  and 

n(w)  oc  a/-1  oc  |fc|-Q , 

(4.3) 
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which  stem  from  the  more  general  Rayleigh-Jeans  distribution 


tcrjM  = 


ci 


C2  +  W 


(4.4) 


They  correspond  respectively  to  equipartition  of  particle  number  N  and  quadratic  energy  E 


N  = 
E  = 


J  n(k)dk  =  J  N(iv)  dtv , 

J  tv{k)n{k)dk  =  J  tvAf(tv)  dev . 


(4.5) 

(4.6) 


The  cases  y  =  0, 1  give  the  non-equilibrium  Kolmogorov-type  solutions,  respectively 

n(tv)  oc  tv  0/  *  +  7  oc  |fc|”2^/3-1+a/3  ,  (4.7) 

and 

n(u/)  oc  a;  oc  |fc|-2^3-1 ,  (4.8) 


which  exhibit  typical  dependence  on  the  parameter  j3  of  the  interaction  coefficient.  The  latter 
solutions  are  more  interesting  since  realistic  sea  spectra  are  of  Kolmogorov-type  by  analogy. 
For  the  case  a  =  1/2  and  /?  =  0,  the  Kolmogorov-type  solutions  are 

n(tv)  oc  w-5/3  a  |fc|~5/6 ,  (4.9) 

n(uj)  oc  u/“2  a  |fc|-1 .  (4.10) 

Both  exponents  satisfy  the  conditions  of  locality  (3.10)-(3.11). 


5  Nature  and  sign  of  the  fluxes 

The  stationary  non-equilibrium  states  are  related  to  fluxes  of  integrals  of  motion,  namely  the 
quantities  N  and  E  in  our  four-wave  interaction  problem.  We  define  the  flux  of  particles  (or 
wave  action)  and  energy  as,  respectively 

=  -f  2^*.  (M) 

P(u)  =  (5.2) 

Here,  Eq.  (4.7),  resp.  Eq.  (4.8),  is  associated  with  constant  flux  Qo-  resp.  Pq,  of  particles, 
resp.  energy.  Let  us  now  mention  a  physical  argument  which  plays  a  crucial  role  in  deciding  the 
realizability  of  the  Kolmogorov-type  spectra.  A  more  detailed  justification  is  provided  below 
in  Section  11  -  see  also  [6]  and  [13].  Suppose  that  pumping  is  performed  at  some  frequencies 
around  w  =  u>f  and  damping  at  u  near  zero  and  a ;  3>  w/.  Weak  turbulence  theory  then  states 
that  the  energy  is  expected  to  flow  from  Uf  to  higher  lj's  (direct  cascade  with  P0  >  0)  while 
the  particles  mainly  head  for  lower  u>’s  (inverse  cascade  with  Qo  <  0).  Accordingly,  we  need  to 
evaluate  the  fluxes  in  order  to  select,  among  the  rich  family  of  power  laws  (4.7)  and  (4.8),  those 
which  are  likely  to  result  from  numerical  simulations  of  Eq.  (2.1)  with  damping  and  forcing 
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By  inserting  Eq.  (3.8)  into  Eqs  (5.1)  and  (5.2),  we  obtain 


oj-v  u~y+1 

Qo  oc  lim - 1 ,  Po  oc  lim  — - 1 , 

2/->0  y  y-+ 1  1  -  y 


which  become 


Qo  oc 


<9/ 

dy 


„  91 

)  Po  OC  — 

»=o  ^ 


y=i 


Using  Eq.  (3.9),  the  derivatives  in  Eq.  (5.4)  can  be  expressed  as 

=  [  £(6.6.6)  (1 +  67  -  6  7  -  67)  61  +  6  -  6  -  6) 

y=0  J  A 

X  <K6“  +  6“  +  6“  ~  1)<^1^2^3) 


a/ 

ay 


a/ 

ay 


s/=l 


[  S (6.6)6)  (1 +  67  -  6 7  -67)  <Ki  +  6  -  6  -  6) 

J  A 

x  [6  In  (1/6)  +  6  In  (1/6)  -  6  In  (1/6)] 
x  5(6“  +  6“  +  6“  -  i)  ^6^6^6 . 


with 


47T 


S  (6.6.6)  =  -7  (666)  “ 


S/2+l_! 


or* 


The  sign  of  each  integral  above  is  determined  by  the  factor  (see  [12]) 

/(  t)  =  i  +  67  -  6 7  -  67  • 

It  is  found  that  /( 7)  is  positive  when 


(5.3) 


(5.4) 


7  <  0  or  7  >  1 .  (5.5) 

For  the  same  values  of  /3  as  those  considered  by  MMT  and  the  additional  value  f3  —  +3,  Table 
1  displays  the  corresponding  frequency  slopes  from  Eqs.  (4.7),  (4.8)  and  the  signs  of  Qo,  Po 
according  to  the  criterion  (5.5). 


-1 

-3/4 

-1/2 

-1/4 

0 

+3 

7 Q 

1/3 

2/3 

1 

4/3 

5/3 

17/3 

sign  of  Qo 

+ 

+ 

0 

- 

- 

7 p 

2/3 

1 

4/3 

5/3 

2 

6 

sign  of  P0 

0 

+ 

+ 

+ 

+ 

Table  1:  Signs  of  the  fluxes  for  the  Kolmogorov-type  solutions. 

Our  calculations  show  that  WT  theory  should  work  most  successfully  for  /3  =  0  (instead 
of  P  =  “I  in  [6])  at  which  they  yield  both  Q0  <  0  and  P0  >  0.  Incidentally,  MMT  reported 
the  smallest  difference  between  numerics  and  theory  for  (3  =  0.  The  cases  with  spectral  slopes 
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less  steep  than  the  Rayleigh- Jeans  distribution  (i.e.  7  <  1)  are  non-physical.  At  best,  a 
thermodynamic  equilibrium  is  expected  in  the  conservative  regime.  Hence,  we  cannot  strictly 
rely  on  the  Kolmogorov-type  exponents  for  /?  =  -1,  -3/4  to  compare  with  the  numerical  results 
in  forced  regimes.  Note  that  for  ft  =  -1/2,  although  we  find  Pq  >  0,  a  pure  thermodynamic 
equilibrium  state  (i.e.  7  =  1)  is  predicted  instead  of  the  inverse  cascade.  This  is  however  not 
valid  because  of  the  necessity  for  a  finite  flux  of  particles  towards  uj  =  0.  The  direct  cascade 
may  then  be  influenced  one  way  or  another,  possibly  making  the  theory  not  applicable  to  the 
whole  spectrum.  Using  both  criteria  (5.5),  we  deduce  that  the  fluxes  of  particles  and  energy 
simultaneously  have  the  correct  signs  in  the  region  of  parameter 


P<~\  and  P>2a~\, 

(5.6) 

3  11 

P  <  ~2  and  P  >  ~2  if  a  =  2  ‘ 

(5.7) 

Since  the  strength  of  nonlinearity  decreases  with  /?,  the  case  j3  <  -3/2,  which  is  close  to  a 
linear  problem,  is  not  that  interesting  from  a  general  viewpoint  and  may  raise  some  difficulties 
in  numerical  studies. 

Restricting  again  to  a  =  1/2  and  j3  =  0,  one  has  for  the  spectrum 


n(u)  =  aP1/3  u)~2 , 


where  P  is  the  flux  of  energy  towards  high  frequencies  and 


a  = 


(5.8) 


(5.9) 


is  the  Kolmogorov  constant.  Numerical  calculations  give  for  a 


a  =  0.376.  (5.10) 

An  important  question  is  the  stability  of  the  stationary  spectra.  This  question  was  studied 
by  Balk  and  Zakharov  in  [14]  from  a  general  point  of  view.  The  particular  situation  discussed  in 
the  present  paper  requires  an  additional  study  based  on  the  work  [14].  However,  one  should  note 
that  instability  of  the  present  spectra  is  unlikely.  The  reason  is  that  the  stationary  spectra  are 
solutions  of  the  kinetic  equation,  which  is  not  sensitive  to  changing  the  sign  of  the  nonlinearity 
in  the  dynamical  equation.  In  other  words,  if  the  Kolmogorov  solution  was  unstable,  it  would  be 
unstable  in  both  cases.  Since  we  observe  the  Kolmogorov  spectrum  in  the  numerical  simulation 
for  one  of  the  signs  of  nonlinearity,  instability  is  unlikely. 

6  Solitons  and  quasisolitons 

Besides  random  radiative  waves,  solitons  are  the  most  interesting  features  of  nonlinear  Hamil¬ 
tonian  models  such  as  the  focusing  NLS.  These  localized  coherent  structures  can  naturally 
emerge  and  persist  as  the  result  of  the  stable  competition  between  nonlinear  and  dispersive 
mechanisms.  It  is  known  that  they  act  as  statistical  attractors  to  which  the  system  relaxes  and 
they  can  influence  the  dynamics  in  a  substantial  way. 
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Equally  important  coherent  structures  are  quasisolitons.  They  could  be  defined  as  solitons 
having  finite  but  long  enough  life  time.  Solitons  and  quasisolitons  can  be  compared  with  stable 
and  unstable  elementary  particles.  Formally,  both  solitons  and  quasisolitons  are  defined  as 
solutions  of  Eq.  (2.3)  of  the  form 

Mt)  =  ei(-n~kV)t  4  •  (6.1) 

Here  and  V  are  constants.  In  the  ar-space, 

ip(x,t)  =  elQt  £(x  -  Vt) ,  (6.2) 

where  £(x)  is  the  inverse  Fourier  transform  of  4>k  and  V  is  the  soliton  velocity.  The  amplitude 
(j>k  satisfies  the  integral  equation 

<j> k  =  J  Tl23k  ^^2^3  +  k2-k3-k)  dkidk2dki  .  (6.3) 

The  “classical”  or  “true”  soliton  is  a  localized  solution  of  Eq.  (6.3).  In  this  case, 

|£(:r)|2 0,  as  |x|->oo.  (6.4) 

This  implies  that  fa  is  a  continuous  function  which  has  no  singularities  for  real  k.  Thus  the 
denominator  in  Eq.  (6.3)  should  not  vanish  on  the  real  axis 

fi  -  kV  +  u;(fc)  ^  0 ,  -oo  <  k  <  -boo  .  (6.5) 

For  w(k)  =  |A;|Q  and  a  <  1,  the  last  condition  is  violated  for  any  V  ^  0.  So  “true  solitons”  can 
exist  only  if  V  =  0. 

Next  we  show  that  “true”  solitons  can  only  exist  for  A  =  -1.  Eq.  (6.3)  can  be  rewritten 
in  the  variational  form 

5{H  +  QN)  =  0.  (6.6) 

Obviously,  fi  >  0  should  hold  (otherwise,  the  denominator  (6.5)  has  zeroes).  Since 

Tusk  =  A  \kik2k3kf^  ,  A  =  ±1 ,  (6.7) 

the  Hamiltonian  is  positive  for  A  =  +1  and  the  condition  (6.6)  can  be  achieved  only  if  ^  =  0. 
There  are  no  solitons  in  this  case.  Meanwhile,  solitons  can  exist  for  A  =  —  1.  A  rigorous  proof 
of  existence  is  beyond  the  frame  of  this  article. 

Quasisolitons  are  a  more  sophisticated  object.  Let  us  allow  the  denominator  (6.5)  to  have 
a  zero  at  k  =  ko  and  suppose  that  ^  is  a  function  which  is  sharply  localized  near  the  wave 
number  k  —  km.  Let  the  width  of  ^  near  k  =  km  be  k.  One  can  introduce 

T{k)  =  j  Tusk  +  &2  ”  ^3  -  k)  dkxdk2dkz  .  (6.8) 

We  might  expect  that 

T(k0)  ^  e-c'-^  \j>(km)\2j>{km) .  (6.9) 

In  other  words,  <j>k  has  a  pole  at  k  =  k-Q  but  the  residue  at  this  pole  is  exponentially  small. 
It  means  that  the  soliton  (6.2)  is  not  exactly  localized  and  goes  to  a  very  small-amplitude 
monochromatic  wave  with  wave  number  k  =  ko  as  x  ->•  — oo. 

If  one  eliminates  the  pole  from  <,6^,  one  gets  a  quasisoliton,  which  is  a  stationary  solution 
of  (2.3)  only  approximately.  Such  a  quasisoliton  lives  for  a  finite  time.  If  this  time  is  long 
enough,  the  quasisoliton  could  become  the  basic  unit  of  wave  turbulence.  This  is  what  we 
believe  may  happen  in  the  MMT  model  with  positive  nonlinearity. 
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7  Soliton  stability  and  collapse 


Coherent  structures  can  play  a  role  in  wave  turbulence  only  if  they  are  stable.  For  A  =  —  1,  a 
soliton  satisfies  the  equation 


(O  +  \k\a)  4>k  =  J  \kik2k3kf/4  fafafc  S{h  +  k2-k3-k)  dhdk2dk3 .  (7.1) 


The  free  parameter  Q  can  be  eliminated  by  the  scaling 


4>k  =  ft  x  (ft  “  &)  , 


(7.2) 


where  x(£)  satisfies  the  equation 


(1  +  m  x(S)  =  J  Itibbtf'4  X1X2X3  m  +6-6-0  dtidhdb .  (7.3) 

Let  us  calculate  the  total  wave  action  in  the  soliton 

N  =  J I4I2  dk  =  CT^r1  No ,  (7.4) 

where 


N0  = 


(7.5) 


The  stability  question  can  be  answered  by  computing  dN/dCl.  As  is  well-known  (see 
[15]),  a  soliton  is  stable  if  dN/dCl  >  0.  In  our  case, 


dN___  (P-a+  1\  N 
dCl  \  a  )  Cl  ' 


(7.6) 


The  soliton  is  stable  if 


P<a-  1 ,  (7.7) 

otherwise  the  soliton  is  unstable.  For  a  =  1/2,  the  condition  of  soliton  instability  reads 

P>~\-  (7.8) 

This  condition  is  satisfied  in  all  the  cases  we  studied. 

The  soliton  instability  leads  us  to  guess  that  the  typical  coherent  structure  in  the  case  of 
negative  nonlinearity  is  a  collapsing  singularity.  Typically,  the  formation  of  such  singularities 
is  described  by  self-similar  solutions  of  the  initial  equations.  Eq.  (2.3)  has  the  following  family 
of  self-similar  solutions 


t)  =  (to  -  t)p+i£  X  [k(to  ~  t)  i ]  ,  (7.9) 

where  p  =  /3~2<^+2  and  e  is  an  arbitrary  constant.  x(£)  satisfies  the  equation 

i  (P  +  ie)  X  +  ^  S  x'  +  ISP  X  +  A  J IS1S2 6SP/4  X1X2X3  *(6  +  6  -  6  -  0  =  0 . 


(7.10) 
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The  soliton  (7.9)  should  stay  finite  when  t  — >  to-  This  requirement  imposes  the  following 
asymptotic  behavior  on  %(£) 

£-»0.  (7.11) 

At  time  t  =  to,  Eq.  (7.9)  turns  to  the  powerlike  function 

-4  CAT" ,  ,  =  ^-^  +  2.  (7.12) 

In  reality,  the  self-similar  solution  is  realized  in  x-space  in  a  finite  domain  of  order  L.  Hence 
the  solution  (7.12)  should  be  cut  off  at  k  ~  1/L.  In  fc-space,  Eq.  (7.9)  represents  the  formation 
of  a  powerlike  “tail”  (7.12).  The  wave  action  concentrated  in  this  tail  must  be  finite.  Therefore 
the  integral 

r  00 

/  \i>k\2dk,  (7.13) 

Jo 

should  converge  as  k  — >  oo.  It  leads  to  the  condition  on  parameters 

(3  >  a  —  1,  (7.14) 

which  coincides  with  the  condition  for  soliton  instability. 

Let  us  plug  (7.9)  into  the  Hamiltonian  in  Fourier  space 

H  =  Ju(k)\i>k\2dk  +  j T\2zk  i>i  foi’-si’t  $(ki  +  fa  -  fa  -  k)  dfadfadfadk , 

=  (to-O^Ho, 

where 

Ho  =  J  16“  Ixl2  d£  +  \J  |6666'?/4  X1X2X3X*  W  1  +6-6-0  d^dfadfadt .  (7.15) 

If  o'  —  1  <  /3  <  2a  —  1,  then  H  — >  00  as  t  -»  to,  unless 

H0  =  0.  (7.16) 

Apparently,  this  condition  can  be  satisfied  only  for  A  =  -1  (negative  nonlinearity).  The  condi¬ 
tion  (7.16)  imposes  implicitly  a  constraint  on  the  constant  e.  In  fact,  it  can  be  realized  only  at 
one  specific  value  of  e,  which  is  an  eigenvalue  of  the  boundary  problem  (7.10)  with  the  boundary 
conditions 

S+2a_2,  £-40, 
x(0  ->  00 ,  16-4  00. 

In  the  case  /?  >  2a  —  1,  H  -4  0  as  t  -4  to.  There  is  no  limitation  on  the  value  of  Hq  and  the 
singularity  can  take  place  for  either  sign  of  A.  If  v  <  1  in  Eq.  (7.15)  or  a  - 1  <  /?  <  a,  a  collapse 
is  the  formation  of  an  integrable  singularity  in  rr-space.  If  v  >  1  or  /3  >  a,  the  singularity  is  the 
formation  of  a  discontinuity  of  the  function  ip(x)  or  its  derivatives. 

The  formation  of  singularities  leads  to  the  formation  in  k- space  of  a  powerlike  spectrum 

nk  a  \ipk\2  ~  | k\~2v  ~  \k\-0+a~2 .  (7.17) 

For  a  =  1/2  and  /?  =  0,  one  obtains 

nk  —  \k\~3/2  ~  u>-3 .  (7.18) 

This  spectrum  can  be  called  Phillips  spectrum  by  analogy  to  the  well-known  “a;-5  spectrum” 
for  deep  water  waves.  As  ui  -4  oo,  it  decays  faster  than  Kolmogorov  spectra. 
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8  More  on  quasisolitons 

Let  us  consider  again  the  case  of  negative  nonlinearity  A  =  —  1  and  denote 

f  =  -n  +  kv  -  w(&)  =  -n  +  kv  -  \k\a .  (8.1) 

If  V  =  0  and  Q,  >  0,  |F|  has  a  minimum  at  k  =  0.  The  Fourier  transform  of  the  solution  is 
concentrated  near  this  minimum  in  a  domain  of  width 

k  ~  fi1/a  .  (8.2) 

Assuming  that  the  soliton  is  smooth  in  x-space,  ^  decays  very  fast  outside  of  the  domain  (8.2). 
So  far  we  assumed  that  V  =  0.  Let  now  V  be  positive  but  very  small.  Then  the  denominator 

i 

F  has  a  zero  at  k  =  ko  —  Va~1.  For  small  V,  the  wavenumber  ko  is  much  larger  than  k  and 
this  zero  occurs  very  far  from  the  domain  which  supports  the  soliton.  This  means  that  fa  has 
a  pole  at  k  =  fco,  but  the  residue  at  this  pole  is  very  small.  The  presence  of  this  pole  means 
that  the  stationary  solution  (6.3)  looks  in  the  x-space  like  a  soliton,  which  is  not  completely 
localized.  As  x  — >  +oo,  it  becomes  a  monochromatic  wave  with  wave  number  ko  and  negligibly 
small  amplitude. 

If  this  “wave  tail”  is  cut  off  in  the  initial  data,  one  has  a  “quasisoliton”  which  slowly 
decays  due  to  radiation  of  energy  in  the  right  direction.  If  V  is  small  enough,  the  lifetime  of 
the  quasisoliton  is  very  long  and  its  shape  is  close  to  the  shape  of  “real”  solitons. 

It  is  unlikely  that  quasisolitons  play  an  important  role  in  wave  turbulence  at  negative 
nonlinearity.  If  V  is  not  small,  their  lifetime  is  too  short;  if  V  is  small,  they  are  unstable  like 
real  solitons.  Quasisolitons  are  more  relevant  in  the  case  of  positive  nonlinearity  A  =  +1. 

Let  us  choose  an  arbitrary  k  =  km  >  0  and  plug  in  Eq.  (6.3) 

v  =  a  ,  tl  =  -(1  -  a)  k°  -  \  a  (1  -  a)  Aft"2  q 2  .  (8.3) 

Then 

F  =  C  -  1*1“  +  a  C1  (*  -  km)  +  \a  (1  -  a)  q 2 .  (8.4) 

Note  that  if  a  <  1  then  F  has  a  zero  at  k  —  ko  <  0  for  any  km.  Hence,  1/F  always  has  a  pole 

on  the  negative  real  axis,  and  the  soliton  (6.3)  cannot  be  a  real  soliton.  But  if  q2  «  Aft,  1/F 

has  a  sharp  minimum  at  k  ~  km.  Introducing 

k  —  \k  kffi | , 

one  has  approximately 

F  a  (1  —  a)  A ft-2  [k2  +  q2]  , 

and  one  gets  for  the  width  of  the  maximum 

Kctq.  (8.7) 

If  k  <  |&o|,  one  can  construct  a  quasisoliton  which  is  supported  in  fc-space  near  km.  In  the 
general  case,  |fc0|  a  km.  lfa  =  1/2  and  q  =  0,  one  can  easily  find 

fco  =  ~(\/2  —  l)2  km .  (8.8) 


(8.5) 

(8.6) 
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The  quasisoliton  moves  to  the  right  direction  with  the  velocity  V(km)  and  radiates  back¬ 
ward  monochromatic  waves  of  wavenumber  ko.  The  shape  of  the  quasisoliton  can  be  found 
explicitly  in  the  limit  q  -¥  0.  Now  k  <  km  and  one  has  approximately 

J  \k\k2kzk\^^  +  k2  —  k3  —  k)  dk\dk2dk3  (8-9) 

—  km  J  <Kk1  +  «2  -  «3  -  k)  dKidKidnz  . 

Taking  into  account  Eq.  (8.6),  one  can  rewrite  Eq.  (6.3)  as 

1  a  (1  -  a)  A£>-2(k2  +  q2)  <j>K=k!}n  J  fafaH  <S(«i  +  «2  -  «3  -  «)  dKidK2dK3  .  (8.10) 

With  the  help  of  inverse  Fourier  transform,  one  can  transform  (8.10)  into  the  stationary  NLS 


-  a  (1  -a)k* 


-2 


dx2 


+  q2(p 


=  ki\<t>\2</>, 


(8.11) 


which  has  the  soliton  solution 


<j)(x)  = 


a  (1  —  a)  q 
&ro~“+2  cosh9x 


It  gives  the  following  approximate  quasisoliton  solution  of  Eq.  (2.1)  with  A  =  +1: 
V>(x,t)  =  4>{x  -  Vi)  eint+ik™(x-vt) , 

n  =  -(l-a)fc“  -\a{l -a)k^~2q2, 

V  =  ak*-1. 


(8.12) 


(8.13) 


The  quasisoliton  (8.13)  is  an  “envelope  soliton”,  which  can  be  obtained  directly  from  Eq.  (2.1). 
Simply  inject 


$(x,t)  =  U(x,t)  ? 


(8.14) 


and  use  the  binomial  expansion 

\k\aU  +  a  |*r 1  (-i£)  U  +  \  a  (o  -  1)  |Af-2  (-i£)  %  U  +  •  •  •  . 

(8.15) 

Plugging  Eq.  (8.15)  into  Eq.  (2.1)  with  A  =  +1,  one  obtains  a  differential  equation  of  infinite 


order 

•( 

lr+yir)  ~  l^+l3u+---, 

(8.16) 

V  =  ak^1. 

Here 

LlV  =  5 

f&TT 

a(l-a)k^2-^  +  ki\U\2U, 

(8.17) 

L3U  =  i 

j  “(«  -  *)(“  -  2>  fcsr3  0  -e'&'wr  % 

,  (8.18) 

LiU  =  •• 

(8.19) 

d 

a 

gikx  jj _  gihx 

dx 
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Taking  into  consideration  only  the  first  nontrivial  term  L2  17,  one  gets  the  nonstationary  NLS 

‘(f +vS)=5a<i-")t»'2S+^|c,|!c'-  (8-2o) 

It  has  a  soliton  solution 

U(x,t)  =  <j>{x  -  Vt )  e-|i«(i-a)fcr2</2t .  (8.21) 

To  find  the  shape  of  the  quasisoliton  more  accurately,  one  should  keep  in  the  right  hand  side  of 
Eq.  (8.15)  a  finite  (but  necessary  odd!)  number  of  terms.  The  expansion  in  Eq.  (8.16)  runs  in 
powers  of  the  parameter  q/km.  Note  that  one  cannot  find  the  lifetime  of  the  quasisoliton.  The 
lifetime  grows  as  e^k°^q  and  its  calculation  is  beyond  the  perturbation  expansion. 

As  a  matter  of  fact,  the  parameter 


is  crucial  for  quasisolitons.  The  smaller  it  is,  the  closer  the  quasisoliton  is  to  a  “real  soliton1’. 
The  amplitude  of  a  quasisoliton  is  proportional  to  e.  Quasisolitons  of  small  amplitude  satisfy 
the  integrable  NLS  and  are  stable.  It  is  not  obvious  for  quasisolitons  of  finite  amplitude.  One 
can  guess  that  at  least  in  the  case  /?  >  0,  when  collapse  is  not  forbidden,  there  is  a  critical 
value  of  the  amplitude  of  a  quasisoliton  ec  such  that  for  e  >  ec  it  is  unstable  and  generates  a 
singularity  at  a  finite  time.  Our  numerical  experiments  confirm  this  conjecture  for  /3  =  +3. 

Quasisolitons  move  with  different  velocities  and  collide.  If  the  amplitudes  of  the  quasisoli¬ 
tons  are  small  and  their  velocities  are  close,  they  obey  the  NLS  and  their  interaction  is  elastic. 
One  can  guess  that  the  same  holds  for  small-amplitude  quasisolitons  even  if  their  velocities  are 
quite  different.  This  is  not  obvious  for  quasisolitons  of  moderate  amplitude.  One  can  think 
that  their  interaction  is  inelastic  and  leads  to  the  merging  and  formation  of  a  quasisoliton  of 
larger  amplitude. 

9  Nonlinear  frequency  shift 

Let  us  consider  one  more  important  nonlinear  effect.  In  a  linear  system,  the  harmonic  of 
wavenumber  k  oscillates  with  the  frequency  In  the  presence  of  nonlinearity,  the  frequency 
changes  due  to  the  interaction  with  other  harmonics.  In  a  weakly  nonlinear  system,  the  fre¬ 
quency  is  modified  by  a  functional  depending  linearly  on  the  spectrum 


w(fc)  -4  w(k)  +  J  T\k  ni  dki . 

(9.1) 

It  is  easy  to  show  that  Ty. 

can  be  expressed  in  terms  of  the  coefficient  in  Eq.  (2.3) 

as 

Tik  =  22ifcifc . 

(9.2) 

For  the  MMT  model, 

Tlk  =  2X(hk)^2. 

(9.3) 

For  0  =  0, 

T\k  =  2X  =  ±2, 

(9.4) 

16 


and 


w±( A; )  =  u}(k)  ±2  N , 


(9.5) 


where  N  =  f  \i>k\2dk  is  the  total  number  of  particles. 

In  the  general  case  0  ±  0,  renormalization  of  the  frequency  leads  to  modified  resonance 
conditions  (2.8)-(2.9).  But  in  the  particular  case  —  0,  renormalization  terms  in  Eq.  (2.9) 
cancel  and  the  resonance  conditions  in  the  first  nonlinear  approximation  remain  unchanged.  At 
the  same  time,  the  difference  of  frequencies  for  different  signs  of  nonlinearity  has  the  form 

u+{k)-u~(k)  =  4N.  (9.6) 

In  our  case,  it  does  not  depend  on  the  wave  number. 

10  On  the  MMT  model  spectrum 

In  [6],  MMT  found  that  in  the  case  of  positive  nonlinearity  the  spectrum  of  wave  turbulence  is 
well  described  by  the  formula  (MMT  spectrum) 


njk  ~  AT^-1 .  (10.1) 

They  checked  this  result  for  a  =  1/2  and  different  values  of  (3.  Our  experiments  are  in  agree¬ 
ment  with  (10.1).  In  [7],  it  was  found  that  the  MMT  spectrum  can  appear  for  either  sign  of 
nonlinearity.  So  far  there  is  no  proper  theoretical  derivation  of  the  MMT  spectrum.  In  this 
section,  we  offer  some  heuristic  derivation  of  (10.1). 

Assuming  formula  (3.2)  to  be  exact,  the  problem  of  closure  for  the  equation  on  particle 
number  lies  in  the  expression  of  Im  Jmk  in  terms  of  n*.  This  expression  should  a  priori  satisfy 
the  conditions  of  symmetry 


Im  Ji23fc  =  Im  J2i3fc  =  Im  Ji2fc3  =  -Im  J3fci2  • 


(10.2) 


Moreover,  one  can  assume  that  the  nonlinearity  is  weak  and  that  the  wave  energy  is  roughly 


(10.3) 


Prom  conservation  of  energy,  one  obtains 


J  Tl23k  (wi  4-  U2  -  a>3  -  w)  Im  J123fc  dkidk2dk3dk  =  0 .  (10.4) 

Hence  one  must  have 


Im  Jl23k  —  <5(^1  +  w2  —  —  u>) .  (10.5) 

For  Gaussian  wave  turbulence,  the  real  part  of  Jl23k  is  given  by  Eq.  (3.3)  and  dimensional 
analysis  gives 


He  J123*; 


k  ‘ 


(10.6) 
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Up  to  this  point,  our  consideration  was  more  or  less  rigorous.  Now  we  present  a  heuristic 
conjecture.  We  suppose  that  the  imaginary  part  of  the  four-wave  correlator  has  the  same 
scaling  as  the  real  part.  In  other  words,  it  is  quadratic  in 

If  one  takes  into  account  the  necessary  conditions  (10.2),  (10.4)  and  the  scaling  (10.6)  for 
Im  Ji23fc?  there  are  not  that  many  possibilities  for  the  construction  of  Im  Just-  We  offer  the 
following  closure 

_  _  f  dcui  doj2  dcj 3  du>\  w  x 

Im  J123A:  =  a  (  0^  +  Qfo  +  +  fa)  S  (wi  +  ^2  “  -  a;)  (n  in2  -  nsnk)  ,  (10.7) 

where  a  <C  1  is  a  dimensionless  constant.  The  closure  leads  to  the  kinetic  equation 

w  =  2o/Tl23fc  (nin2  - n3Hk)  (l!r +  +  !§  +  SO  (10-8) 

x  5  (a;i  +  0J2  —  c^3  “  w)  5  {k\  +  /c2  —  k%  -  k)  dk\dk2dk$  . 

It  is  easy  to  check  that  the  Kolmogorov  solution  of  Eq.  (10.8)  leads  to  the  MMT  spectrum. 
Eq.  (10.8)  resembles  the  Boltzmann’s  equation  for  interacting  particles.  Apparently,  it  can  make 
sense  only  if  aT\23fc  >  0.  Otherwise,  the  iif-theorem  and  the  second  law  of  thermodynamics  will 
be  violated.  We  must  stress  that  the  formula  (10.7)  is  heuristic  and  has  no  rigorous  justification. 


11  Particle  and  energy  balance 


In  the  presence  of  damping  and  linear  instability,  Eq.  (2.3)  can  be  written  in  the  form 


where 


1  dt 


777  +i  D(k)i>k, 

srk 


H  =  J  u{k)  h h?dk+l-  J  Ti23k  V’lV’aV’alV’fc  £(fci  +  k2-fa-  k)  dkidkidk^dk , 


(11.1) 


(11.2) 


is  the  Hamiltonian,  D(k)  is  the  damping  or  the  growth  rate  of  instability  depending  on  its  sign. 

Let  N  =  f  |^Sfc|2d&  be  the  total  number  of  particles  in  the  system.  From  (11.1),  one  can 
obtain  the  exact  equation  for  the  particle  balance 


^  =  Q  =  2  j  D{k)m2dk.  (11.3) 

After  averaging,  one  has 

®  =  2  j  D{k)nkdk  =  {Q).  (11.4) 

The  flux  of  particles  Q  is  a  linear  functional  of  nk  at  any  level  of  nonlinearity. 

For  the  energy  flux,  one  has  the  exact  identity 

dH  r 

P  =  —  =  2J  u(k)D(k)  \ipk\2dk 

+  2  J  +  P{k2)  +  D{ka)  +  D(k)]  T\23k  V'lltaV’aiV’fc 

x  6{k\  +  k2  —  k$  —  k)  dkidk2dk^dk .  (11.5) 
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For  the  averaged  density  of  energy,  one  has 

(P)  =  2  J  u{k)D(k)nkdk 

+  \  I  [D(k i)  +  D(k2)  +  D(k3)  +  D(k)]Tmk  Re  Jmk 
x  6(ki  +  &2  —  —  k)  dk\dk2dk^dk  .  (11.6) 

Assuming  that  Gaussian  statistics  holds,  one  can  write 

Re  J123A:  -  nin2  [5{h  -  k3)  +  S(ki  -  k )]  ,  (11.7) 

and  one  obtains  after  simple  calculations 

(P)  =  2  J  oi(k)  D(k)nkdk ,  (11.8) 

where  L o(k)  =  u j(k)  +  f  T\k  n\  dk\  is  the  renormalized  frequency. 

In  the  case  /?  =  0  and  =  ±2, 

(P)  =  2  J  w(Jfe)  D(k)  nk  dk  +  2  A  N  (Q) .  (11.9) 

In  the  stationary  state,  ( Q )  =  0,  ( P )  =  0  and  the  balance  equations  are 

J  D{k)rikdk  =  0,  (11.10) 

J  w(k)  D(k)  rik  dk  —  0.  (11.11) 

In  this  particular  case,  renormalization  of  the  frequency  does  not  influence  the  balance  equa¬ 
tions. 

The  balance  equations  (11. 10)- (11. 11)  can  be  rewritten  as 

Qo  =  Q+  +  Q“,  (11.12) 

Po  -  P+  +  P”,  (11.13) 

where  Po  and  Qo  are  the  input  of  particles  and  energy  in  the  area  of  instability  cv  ~  u/o-  Q+ 
and  P+  are  the  sinks  of  particles  and  energy  in  the  high  frequency  region  lj  ~  <*;+.  Q~  and  P” 
are  the  sinks  in  the  low  frequency  region  u  ~  cj~  . 

Roughly  speaking, 


Po 

-  ^oQo, 

(11.14) 

p+ 

^  w+Q+, 

(11.15) 

p- 

-  v-Q-, 

(11.16) 

and  the  balance  equations  can  be  written 

as 

Qo  = 

Q+  +Q~ , 

(11.17) 

w0Qo  =*  w+Q+  +  u  Q  .  (11.18) 
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Hence 


Q+  wo  -  w  P+  w+  wo  -  w 

<2^  _  _  _  _ 

Q~  u)+  -  ujq  ’  P~  a;-  a;+  -  cj0  * 


For  a;  ~  u>o  ^  or1-,  one  has 


Q+  ^  u>o  -  a;' 
~  o;+ 


P+  a?o  —  w 
____  _ 

P-  W 


(11.19) 


(11.20) 


In  other  words,  if  wo  <C  w+,  almost  all  particles  are  absorbed  at  low  frequencies.  The  amounts 
of  energy  absorbed  in  both  ranges  have  the  same  order  of  magnitude.  These  conclusions  are 
valid  only  under  the  hypothesis  of  approximate  Gaussianity  of  wave  turbulence. 
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and  for  the  quadratic  part  of  energy 


N 

2 

E=  v{k„)  \i>n\2 .  (12.5) 

n=-f+l 

The  linear  frequency  term  is  treated  exactly  by  an  integrating  factor  technique,  removing 
it  from  the  timestepping  procedure.  As  emphasized  by  MMT,  we  thus  avoid  the  natural  stiffness 
of  the  problem  as  well  as  possible  numerical  instabilities.  Consequently,  we  do  not  need  to 
shorten  the  inertial  interval  by  downshifting  the  cutoff  of  ultraviolet  absorption  (as  in  [4]).  The 
nonlinear  term  is  calculated  through  the  Fast  Fourier  Transform  by  first  transforming  to  real 
space  where  a  multiplication  is  computed  and  then  transforming  back  to  spectral  space.  For 
the  multiplication  operation,  twice  the  effective  number  of  grid  points  are  required  in  order  to 
avoid  aliasing  errors.  A  fourth-order  Runge-Kutta  scheme  integrates  the  conservative  model  in 
time,  giving  a  solution  to  which  the  diagonal  factor 

e[F(k)+D(k)]At 


is  applied  at  each  time  step  At. 

13  Numerical  results  for  (3  =  0,  A  =  ±1 

A  series  of  numerical  simulations  of  Eq.  (12.1)  with  resolution  up  to  2048  de-aliased  modes  has 
been  performed.  We  choose  the  case  /3  =  0  as  the  candidate  for  testing  weak  turbulence  in  our 
experiments.  Both  cases  A  =  ±1  are  examined,  providing  an  additional  test  of  the  theory,  and 
the  study  is  focused  on  the  direct  cascade.  Forcing  is  located  at  large  scales  and  the  inertial 
interval  is  defined  by  the  right  transparency  window  kf  <  k  <  (where  kf  and  kd  are  the 
characteristic  wave  numbers  of  forcing  and  ultraviolet  damping  respectively).  As  displayed  in 
Table  1,  the  theoretical  spectrum  which  can  be  realized  in  this  window  is 

nk(xk~l.  (13.1) 

Typically,  initial  conditions  are  given  by  the  random  noise  in  the  spectral  space.  Simulations 
are  run  until  a  quasi-steady  regime  is  established  which  is  characterized  by  small  fluctuations  of 
the  energy  and  the  number  of  particles  around  some  mean  value.  Then  time  averaging  begins 
and  continues  for  a  length  of  time  which  significantly  exceeds  the  characteristic  time  scale  of 
the  slowest  harmonic  from  the  inertial  range  (free  of  the  source  and  the  sink).  In  turn,  the 
time- step  of  the  integration  has  to  provide,  at  least,  accurate  enough  resolution  of  the  fastest 
harmonic  in  the  system.  As  our  experiments  show,  one  has  to  use  an  even  smaller  time-step 
than  defined  by  the  last  condition  :  the  presence  of  fast  nonlinear  events  in  the  system  requires 
the  use  of  a  time-step  At  =  0.005,  which  is  40  times  smaller  than  the  smallest  linear  frequency 
period.  Time  averaging  with  such  a  small  time  step  leads  to  a  computationally  time-consuming 
procedure  despite  the  one-dimensionality  of  the  problem. 

From  now  on,  we  will  present  numerical  results  in  the  specific  situations  v~~  —  196.61  (A  = 
±1),  =  5.39  x  10  48  (A  =  H-l)  or  =  2.16  x  10-47  (A  =  —1),  and  fj  =  0.2,  nonzero  only 

for  kj  e  [6, 9]  (A  =  ±1). 

The  numerical  simulations  clearly  display  the  development  of  dynamical  chaos  and  sta¬ 
tistically  uniform  turbulence.  Both  the  amplitude  and  the  phase  of  each  harmonic  fluctuate 
independently  of  each  other.  Fig.  1-4  show  the  behavior  of  the  seventh  and  eighth  harmonics. 
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Fig.  5-8  show  the  behavior  of  the  real  and  imaginary  parts  of  the  amplitude  of  the 
harmonic  k  =  200.  One  sees  amplitude-modulated  oscillations  with  carrying  frequency  close  to 
the  corresponding  linear  frequency  of  the  harmonic  w  ~  14. 

Fig.  9-12  represent  Fourier  transforms  in  time  of  the  evolution  of  the  harmonic  k  =  200 
from  the  previous  pictures.  One  can  see  that  the  maximum  of  the  spectra  corresponds  to  the 
linear  frequency  shifted  in  accordance  with  the  nonlinearity  sign  A  =  ±1. 

Fig.  13-14  demonstrate  the  behavior  of  the  fourth  and  sixth-order  moments  as  functions 
of  the  second-order  moment.  They  fit  the  Gaussian  laws  very  well.  They  provide  a  justification 
of  the  initial  conjecture  that  the  statistics  of  the  turbulence  is  close  to  Gaussian. 

Fig.  15  represents  the  time  evolution  of  the  quadratic  energy  E  for  A  =  ±1  with  the 
same  amplitude  of  forcing.  The  curves  are  plotted  over  the  interval  t  E  [5000, 10000]  where 
the  time  averaging  actually  takes  place.  One  obviously  sees  that  the  systems  have  already 
reached  the  steady  state.  Their  energies  moderately  fluctuate  about  mean  values  which  are 
Em  —  19  (A  =  +1)  and  Em  ~  9  (A  =  -1).  This  significant  difference  with  respect  to  the  sign 
of  A  is  quite  unexpected  from  the  viewpoint  of  the  WT  theory  since  the  same  rate  of  forcing  is 
imposed  in  both  systems.  We  can  make  the  same  remarks  about  the  evolution  of  the  number 
of  particles  N.  In  Fig.  16,  the  mean  values  stay  near  Nm  ~  3  (A  =  +1)  and  Nm  ~  1  (A  =  —  1) 
so  that  their  relative  difference  is  even  bigger  than  for  E.  Fluctuations  also  spread  much  more 
in  the  case  A  =  +1. 

In  Fig.  17,  the  stationarity  as  well  as  the  gap  between  both  signs  of  A  are  verified  again 
in  the  time  evolution  of  the  average  nonlinearity  e.  We  define  the  average  nonlinearity  in  the 
system  as  the  ratio  of  the  nonlinear  part  to  the  linear  part  of  the  Hamiltonian  e  =  I-Hatl/UlI, 
each  part  being  calculated  over  the  whole  field.  Of  course,  this  definition  does  not  really 
make  sense  when  external  forces  are  applied  but  it  provides  a  relatively  good  estimation  of  the 
level  of  nonlinearity  once  the  systems  reach  the  steady  state.  Note  here  that  the  mean  values 
em  ~  0.4  (A  =  +1)  and  em  ~  0.2  (A  =  —1)  are  relatively  small.  Thus,  the  condition  of  small 
nonlinearity  required  by  the  theory  holds  for  both  systems.  This  conclusion  is  also  supported 
by  comparing  Fig.  10  and  Fig.  12.  It  is  seen  that  the  difference  of  frequencies  caused  by 
nonlinearity  is  relatively  small.  We  point  out  that  in  our  numerical  experiments  e  could  not  be 
taken  too  small  (that  is,  e  <  10-3)  for  two  reasons.  First,  the  nonlinear  turnover  time  grows 
longer  and  the  energy  flux  is  too  weak  to  act  effectively.  Second,  one  may  catch  the  undesirable 
frozen  turbulence  [10]  due  to  the  disappearance  of  quasiresonances.  One  should  note  that,  in 
general,  frozen  turbulence  arises  more  easily  in  one-dimensional  problems  due  to  fewer  degrees 
of  freedom  than  in  higher-dimensional  problems. 

The  difference  between  the  cases  A  =  ±1  is  especially  conspicuous  if  one  considers  the 
dissipation  rates  of  particles  and  quadratic  energy  in  the  left  transparency  window 

Q~  =  2  v~\k\~d~  \i>k\2dk  ,  P~  =  2  I*' v-\k\~d~  a>(*)|&|2dfc, 

Jo  Jo 

and  in  the  right  transparency  window 

fkmax  ,  A  rkmax  . 

Q+  =  2/  v+\k\d  \A\2dk  ,  P+  =  2/  v+\k\d+  u(k)  \ij>k\2  dk , 

Jkf  Jkf 

where  kmax  corresponds  to  the  highest  mode  in  the  spectrum.  Fig.  18-21  represent  the  time 
evolution  of  these  quantities  and  their  time-averaged  values  are  collected  in  Table  2. 

One  can  see  that  the  case  A  =  +1  quantitatively  fits  WT  theory.  Indeed,  in  this  case 
Q+/Q  —  0.046  1  and  P+/P~  ~  0.94.  But  in  the  case  of  negative  nonlinearity  A  =  —  1  the 
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A 

□ 

E 

Q~ 

Q+ 

P~ 

P+ 

+1 

3 

fil 

0.276 

0.258 

-1 

1 

9 

0.0098 

lilfgffil 

0.014 

1.430 

Table  2:  Time-averaged  values  of  the  wave  action,  quadratic  energy  and  corresponding  fluxes 
in  the  stationary  state. 


situation  is  opposite.  In  this  case  Q^/Q~  —  4.9  and  P+/P”  ~  102  which  means  that  most  of 
both  quadratic  energy  and  particles  are  transported  to  high  frequencies. 

Comparison  of  the  turbulence  levels  and  fluxes  of  particles  Q+  for  both  signs  of  nonlin¬ 
earity  leads  to  a  paradoxal  result.  At  A  =  -1  the  total  number  of  particles  is  three  times  less 
than  at  A  =  +1,  while  the  dissipation  rate  of  particles  is  higher  by  one  order  of  magnitude. 
It  can  be  explained  only  by  the  presence  in  this  case  of  a  much  more  powerful  mechanism  of 
nonlinear  interactions,  which  provides  very  fast  wave  particles  transport  to  high  frequencies. 
In  our  opinion,  this  mechanism  is  wave  collapse,  studied  theoretically  in  Section  7.  Sporadic 
collapsing  events  developing  on  top  of  the  WT  background  could  send  most  of  particles  to  high 
wavenumbers  without  violation  of  energy  conservation,  because  in  each  self-similar  collapse 
structure  the  amount  of  total  energy  is  zero. 

We  observed  such  collapsing  events  in  our  numerical  experiments.  Fig.  22  displays  the 
collapse  event  taking  place  at  the  point  x  =  1.006  at  time  t  =  5000.19.  One  can  conjecture  that 
the  collapses  are  described  by  self-similar  solutions.  For  such  solutions  H  =  0.  It  means  that 
the  collapse  can  carry  particles  to  high  frequencies,  without  carrying  any  energy  at  that  time\ 
As  far  as  the  Hamiltonian  is  the  difference  of  quadratic  and  quartic  terms  and  both  of  them  go 
to  infinity,  it  becomes  possible  to  explain  the  apparent  contradictions  of  the  dissipation  rates. 

The  hypothesis  related  to  the  prevailing  role  of  collapses  at  A  =  —1  is  corroborated  by 
the  following  facts: 

1.  Intermittency  in  dissipation  rates  of  quadratic  energy  and  particles  for  A  =  -1  is  much  higher 
than  for  A  =  +1  in  the  region  of  large  wave  numbers.  This  intermittency  can  be  explained  by 
outbursts  of  dissipation  when  wave  collapses  occur. 

2.  The  analysis  of  time  Fourier  transforms  of  separate  harmonics  (we  take  k  =  200)  shows  the 
presence  of  two  components,  see  Fig.  9.  The  peak  at  u  ~  13  corresponds  to  a  linear  wave  with 
a  moderate  nonlinear  shift  of  frequency.  This  is  the  “weak  turbulence”  component  of  the  wave 
field.  Another  component  is  roughly  symmetrical  with  respect  to  the  reflection  u)  — >  —co  with 
the  maximum  slightly  below  oj  =  0.  This  is  certainly  a  strongly  nonlinear  component  which 
could  be  associated  with  wave  collapses. 

Another  indication  of  the  difference  of  the  wave  dynamics  in  the  cases  A  =  +1  and  A  =  -1 
follows  from  the  following  experiment.  Fig.  23-24  show  the  early  stages  in  the  conservative 
evolution  of  the  same  isolated  initial  condition 

(x— ff)2 

ip(x)  =  ip oe  ,  <7  =  0.5. 

In  the  case  A  =  -1,  a  sufficiently  large  initial  condition  collapses  into  a  sharp  spike,  while  in  the 
case  A  =  +1  it  decays.  This  experiment  could  serve  as  an  evidence  of  the  finite-time  singularity 
formation  for  the  case  A  =  —  1. 

Now  we  discuss  the  stationary  isotropic  spectra  of  turbulence  which  are  displayed  in  Fig. 
25-28.  We  plotted  on  the  same  pictures  the  Kolmogorov  spectra  calculated  by  putting  either 
P  =  P+  =  1-430  (A  =  -1)  or  P  =  P+  =  0.258  (A  =  +1)  and  a  =  0.376  in  Eq.  (5.8).  In  Fig. 
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27-28,  one  can  see  that  for  both  cases  this  spectrum  provides  a  higher  level  of  turbulence  than 
the  observed  one.  For  A  =  -1  this  difference  is  almost  of  one  order  of  magnitude.  For  A  —  +1, 
the  observed  spectrum  almost  coincides  with  the  weak  turbulence  one  at  low  frequencies  and 
then  decays  faster  at  higher  wave  numbers  (approximately  as  MMT  spectrum  in  Fig.  26). 

It  is  interesting  that  for  A  =  -1  the  high  frequency  asymptotics  is  fairly  close  to  the 
one  predicted  by  WT  theory  (Fig.  25).  One  can  explain  this  fact  as  follows.  In  this  case,  the 
turbulence  is  the  coexistence  of  collapsing  events  and  weak  turbulence.  Collapses  carry  most 
of  the  fluxes  of  particles  and  quadratic  energy  to  high  frequencies.  But  their  contribution  to 
the  high-frequency  part  of  the  spectrum  is  weak,  because  they  produce  Phillips-type  spectra, 
decaying  very  fast  as  k  — ^  oo.  In  our  case,  this  spectrum  is 

t ik  ~  AT3/2  .  (13.2) 

Hence  as  k  ->  oo,  only  the  WT  component  survives.  Even  P  ~  10-2  P+  is  enough  to  provide 
an  observable  tail  in  the  WT  Kolmogorov  spectrum. 

We  should  stress  out  again  that  at  A  =  +1  the  picture  of  turbulence  matches  the  WT 
prediction  both  quantitatively  and  qualitatively.  Meanwhile,  the  spectrum  at  high  k's  is  steeper 
and  closer  to  the  MMT  formula.  So  far  we  cannot  give  a  consistent  explanation  of  this  fact. 
We  can  just  guess  that  it  is  somehow  connected  with  quasisolitons.  As  an  illustration,  Fig.  29 
shows  the  conservative  evolution  of  the  initial  quasisoliton  (8.13)  with  parameter  q/km  =  0.1, 
which  is  small  enough  to  justify  the  Taylor  expansion  used  in  its  derivation.  As  expected,  we 
observe  that  the  solution  propagates  and  persists  over  a  relatively  long  time.  This  similarity 
between  quasisolitons  and  real  solitons  is  verified  even  better  in  Fig.  30-33  where  two  initial 
quasisolitons  with  q/km  =  0.2  for  the  smaller  one  and  qjkm  =  0.25  for  the  bigger  one  collide 
almost  elastically.  Note  here  that  the  solution  with  smaller  amplitude  moves  with  a  greater 
velocity. 

14  Numerical  results  for  (3  =  +3  and  A  =  +1 

Another  series  of  experiments  has  been  performed  for  the  case  (3  —  +3  and  A  =  + 1.  This  case 
is  especially  attractive  due  to  the  fact  that  the  intensity  of  interaction  grows  with  characteristic 
wavenumber  in  Fourier  space  and  one  can  expect  reduced  “frozen”  turbulence  effects  compared 
to  the  case  /?  =  0.  Another  motivation  is  the  fact  that  the  scaling  of  the  interaction  kernel 
reproduces  the  kernel  for  gravity  water  waves.  Therefore,  Eq.  (12.1)  with  o:  =  1/2,  /?  =  +3  can 
be  considered  as  a  model  of  turbulence  of  the  ocean  surface. 

The  numerical  simulation  of  Eq.  (12.1)  was  performed  on  a  grid  of  2048  points  in  the 
real  space  domain  of  length  27t.  Parameters  of  the  forcing  are  defined  by 

f  0.001  if  30  <  A:  <42, 

\  0  otherwise , 

and  parameters  of  damping  in  the  “hyperviscosity”  form  by 

{-0.05(fc  —  4)8  if0<fc<4, 

-0.1(fc  -  824)2  if  824  <  A;  <  1024 , 

0  otherwise . 

Aliasing  effects  were  not  of  concern  due  to  the  run-time  control  of  the  fastness  of  the  spectrum 
decay  toward  high  wavenumbers. 
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The  time-step  of  integration  was  equal  to  ^  of  the  inverse  fastest  linear  frequency  in 
the  problem.  Such  a  small  value  was  chosen  due  to  the  fact  that  the  time  dependence  of  the 
individual  Fourier  harmonics  corresponding  to  intermediate  range  wavenumbers  showed  the 
presence  of  processes  of  time  scale  smaller  than  the  smallest  linear  time  in  the  system.  This 
observation  was  an  initial  indication  of  the  significant  role  of  nonlinearity  in  the  problem  under 
consideration. 

Equation  (12.1)  was  integrated  numerically  over  long  times  for  different  kinds  of  initial 
conditions:  low  level  random  noise  and  single  harmonic  excitation  ( k  =  30)  initial  conditions. 
While  initial  stages  of  computations  were  quantitatively  different,  the  later  stages  of  evolution 
were  strikingly  similar.  Starting  from  big  enough  times,  the  wave  system  was  separated  into 
several  soliton-like  moving  structures  and  low-amplitude  quasi-linear  waves.  Processes  of  in¬ 
teraction  of  solitons  and  waves  slowly  redistributed  the  number  of  waves  in  a  way  leading  to 
the  growth  of  initially  bigger  solitons  and  the  collapse  of  initially  smaller  solitons.  Finally  the 
system  was  clearly  separated  into  a  state  with  one  moving  soliton  and  quasi-linear  waves. 

We  interpret  the  observed  phenomenon  as  similar  to  the  “droplet”  effect  observed  earlier 
in  non-integrable  NLS  equation  [16].  The  soliton  solution  turns  out  to  be  the  statistical  attractor 
for  nonlinear  non-integrable  wave  systems:  long  time  evolution  leads  to  the  condensation  of  the 
integral  of  total  number  of  waves  into  the  single  soliton  which  minimizes  the  Hamiltonian. 

Fig.  34-35  show  snapshots  of  the  final  state  of  the  system:  the  single  soliton  is  moving 
with  constant  speed  on  the  background  of  quasi-linear  waves.  A  quantitative  comparison  shows 
that  the  parameters  of  the  observed  object  are  close  to  the  parameters  of  the  quasisoliton 
solution  (8.13). 

One  should  emphasize  that  there  is  a  difference  between  the  situation  observed  in  the 
present  work  and  former  observations  of  “droplet”  effects  in  non-integrable  NLS  equations. 
Solitons  observed  in  [16]  were  exact  stable  solutions  of  the  corresponding  NLS  equation.  Solitary 
solutions  observed  in  the  present  work  are  “quasisolitons”  which  are  unstable  at  least  in  a  certain 
range  of  parameters. 

In  Fig.  36  the  initial  condition  is  the  quasisoliton  (8.13)  with  parameter  q/km  =  0.1.  Here 
again,  it  behaves  as  the  soliton  should:  it  moves  without  any  detectable  change  of  shape.  Fig. 
37  shows  the  evolution  for  qjkm  =  0.3.  One  can  interpret  such  initial  condition  as  a  “deformed” 
quasisoliton.  This  initial  condition  rapidly  develops  moving  singularity  collapsing,  presumably, 
in  finite  time. 

15  Conclusion 

The  MMT  model  with  a  <  1  and  either  sign  of  nonlinearity  exhibits  coherent  structures.  In 
the  case  of  negative  nonlinearity  these  structures  are  weak  collapses.  These  collapses  are  a 
powerful  mechanism  of  energy  dissipation,  which  dominates  in  all  our  numerical  experiments. 
Weak  turbulence  coexists  with  collapses,  and  is  responsible  for  the  formation  of  Kolmogorov- 
type  tails  of  wave  spectra.  But  it  carries  to  high  wave  numbers  just  a  small  part  of  the  energy 
(less  than  5%). 

One  may  hope  to  get  “pure”  weak  turbulence  by  decreasing  the  level  of  nonlinearity.  But 
to  achieve  an  adequate  modeling  of  the  continuous  medium,  one  should  take  a  very  fine  mesh  (at 
least  104  harmonics)  and  apply  forcing  in  a  broad  range  (say  10  <  k  <  100).  Otherwise  effects 
of  “frozen  turbulence”  will  blur  the  picture.  Such  experiments  would  be  very  time-consuming. 

The  case  of  positive  nonlinearity  is  less  clear.  In  this  case  the  picture  of  turbulence  is 
qualitatively  similar  to  weak  turbulence,  but  the  slope  of  the  spectrum  fits  better  the  MMT 


26 


spectrum.  So  far  we  do  not  have  a  satisfactory  explanation  of  this  phenomenon.  Probably  it 
could  be  explained  by  the  presence  of  interacting  quasisolitons.  In  this  case  again,  experiments 
with  a  larger  number  of  harmonics  could  give  a  result  closer  to  WT  predictions. 

The  relative  “suppression”  of  weak  turbulence  in  the  MMT  model  can  be  explained  by 
a  peculiarity  of  the  resonant  conditions.  In  the  one-dimensional  case  with  a  =  1/2,  only  well- 
separated  waves  interact.  Indeed,  one  can  see  from  (2.10)  that 


fa 

fa 


£>0, 


(15.1) 


and  therefore  min  |fe2/fci|  =  9  is  reached  at  £  =  1.  This  phenomenon  can  be  called  “sparsity  of 
resonances” .  Due  to  this  sparsity,  four  wave  resonances  easily  lose  the  competition  with  coherent 
structures  -  collapses  and  quasisolitons.  In  this  sense  the  MMT  model  is  not  an  optimal  object 
for  checking  the  validity  of  WT  theory.  We  can  offer  the  following  model,  which  includes  the 
interaction  of  two  types  of  waves 


H  = 


x 


I  \k\a  ( \ak\2  +  s\bk\2)dk 

j  \kfafafaf^  (<2^a*ci2a3  +  2pi  0^6102^3  +P2  bfrb^fabs) 
S(k  +  fa  -  fa  ~  fa)  dkdfadfadfa  . 


(15.2) 


If  a  >  1  and  /?  <  2a  —  1,  the  corresponding  dynamical  system  does  not  describe  any 
coherent  structures  which  could  compete  with  four- wave  resonances.  Meanwhile,  for  s  ^  1,  it 
describes  nontrivial  resonant  interactions  for  different  waves  propagating  in  the  same  direction. 
The  system  (15.2)  looks  like  a  possible  object  for  the  simulation  of  wave  turbulence.  In  the 
special  case  a  =  2  and  /?  =  0,  it  describes  coupled  NLS  equations. 
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Figure  1:  /?  =  0,  A  =  — 1.  Amplitude  of  the  mode  k  =  7  vs.  time. 


Figure  2:  /3  =  0,  A  =  —  1.  Amplitude  of  the  mode  k  =  8  vs.  time. 
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Figure  5:  =  0,  A  =  —  1.  Time  evolution  of  the  real  part  of  the  amplitude  for  the  mode  k  =  200 

(time  resolution  r  =  0.015). 


Figure  6:  /?  =  0,  A  —  +1.  Time  evolution  of  the  real  part  of  the  amplitude  for  the  mode  k  =  200 
(time  resolution  r  =  0.015). 
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Figure  7:  ft  =  0,  A  =  — 1.  Time  evolution  of  the  imaginary  part  of  the  amplitude  for  the  mode 
k  =  200  (time  resolution  r  =  0.015). 


Figure  8:  /3  —  0,  A  —  +1.  Time  evolution  of  the  imaginary  part  of  the  amplitude  for  the  mode 
k  =  200  (time  resolution  r  =  0.015). 
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Figure  9:  /3  =  0,  A  =  — 1.  Square  amplitude  of  the  Fourier  transform  for  the  mode  k  =  200  vs. 
frequency  (time  resolution  r  =  0.015). 


Figure  10:  /3  —  0,  A  =  — 1.  Same  as  before  but  with  a  zoom  on  a  smaller  frequency  window. 


Figure  13:  /3  =  0,  A  =  +1.  Fourth  (crosses)  and  sixth-order  (circles)  moments  as  functions  of 
the  second-order  moments.  The  straight  lines  are  the  fitted  Gaussian  laws. 


Figure  14:  /?  —  0,  A  =  —  1.  Fourth  (crosses)  and  sixth-order  (circles)  moments  as  functions  of 
the  second-order  moments.  The  straight  lines  are  the  fitted  Gaussian  laws. 


Figure  19:  /3  —  0,  A  =  —1  (left),  A  =  +1  (right).  Dissipation  rate  of  particles  at  high  wave 
numbers  vs.  time. 


Figure  20:  /?  =  0,  A  =  —  1  (left),  A  =  +1  (right).  Dissipation  rate  of  energy  at  low  wave  numbers 
vs.  time. 


Figure  21:  p  =  0,  A  =  -1  (left),  A  =  +1  (right).  Dissipation  rate  of  energy  at  high  wave 
numbers  vs.  time. 
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Figure  23:  /3  =  0,  A  =  — 1.  Evolution  towards  a  collapsing  peak  of  the  isolated  solution  for  the 
initial  amplitude  ^ o  =  2.  Dotted  line  t  =  0,  dashed  line  t  —  0.55,  solid  line  t  =  1.1. 
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Figure  24:  ft  —  0,A  —  +1.  Evolution  towards  decay  of  the  isolated  solution  for  the  initial 
amplitude  =  2.  Dotted  line  t  =  0,  dashed  line  t  =  1.65,  solid  line  t  =  3.85. 


Figure  25:  /3  =  0,  A  —  —  1.  Computed  spectrum  vs.  wave  number.  The  theoretical  slopes  are 
shown  as  well  (A;-1  for  WT  and  fc~5/4  for  MMT). 
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Figure  26:  /?  —  0,  A  =  +1.  Computed  spectrum  vs.  wave  number.  The  theoretical  slopes  are 
shown  as  well  (A;-1  for  WT  and  A;-5/4  for  MMT). 


Figure  28:  /3  —  0,  A  —  +1.  Computed  spectrum  and  WT  spectrum  vs.  wave  number.  The  WT 
spectrum  (straight  line)  is  given  by  n(k)  =  ck~l  with  c  =  a  P1/3  ~  0.24. 


Figure  30:  /3  =  0,  X  —  +1.  Interaction  of  two  initial  quasisolitons  at  t  =  0.  The  smaller  and 
bigger  ones  correspond  to  q/km  —  0.2  and  0.25  respectively. 


Figure  31:  ft  —  0,  A  =  +1.  Interaction  of  two  initial  quasisolitons  at  t  —  37.5. 
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Figure  34:  /?  =  3,  A  =  +1.  Single  moving  soliton,  t  =  61 
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Figure  35:  /3  =  3,  A  =  +1.  Single  moving  soliton,  t  =  10 
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Figure  36:  j3  =  3,  A  =  +1.  Evolution  of  the  initial  quasisoliton  for  q/km  —  0.1.  Solid  line  t  =  0, 
dotted  line  t  —  23.6,  dashed  line  t  —  47.1. 


Figure  37:  /?  =  3,  A  =  +1.  Evolution  of  the  initial  quasisoliton  for  q/km  —  0.3.  Solid  line  t  =  0 
dotted  line  t  =  23.6,  dashed  line  t  =  47.1. 


